Convergence Analysis of Particle Swarm Optimizer and Its Improved Algorithm Based on Velocity Differential Evolution

Abstract

This paper presents an analysis of the relationship of particle velocity and convergence of the particle swarm optimization. Its premature convergence is due to the decrease of particle velocity in search space that leads to a total implosion and ultimately fitness stagnation of the swarm. An improved algorithm which introduces a velocity differential evolution (DE) strategy for the hierarchical particle swarm optimization (H-PSO) is proposed to improve its performance. The DE is employed to regulate the particle velocity rather than the traditional particle position in case that the optimal result has not improved after several iterations. The benchmark functions will be illustrated to demonstrate the effectiveness of the proposed method.

1. Introduction

Algorithms to tackle optimization problems include not only classical techniques such as dynamic programming, branch-and-bound, and gradient-based methods, but also more recent techniques such as metaheuristics [1]. Among the existing metaheuristic algorithms, the particle swarm optimization (PSO) algorithm is a population-based optimization technique developed by Kennedy and Eberhart in 1995 [2]. The PSO has resulted in a large number of variants of the standard PSO. Some variants are designed to deal with specific applications [3–6], and others are generalized for numerical optimization [7–10]. A hierarchical version of PSO (H-PSO) has been proposed by Janson and Middendorf [10]. In H-PSO, all particles are arranged in a tree that forms the hierarchy. A particle is influenced by its own best position and the best position particle in its neighborhood. It was shown that H-PSO performed very well compared to the standard PSO on unimodal and multimodal test functions [10, 11]. H-PSO presents the advantage of being conceptually very simple and requiring low computation time. However, the main disadvantage of H-PSO is the risk of a premature search convergence, especially in complex multiple peak search problems.

A number of algorithms combined various algorithmic components, often originating from algorithms of other research areas on optimization. These approaches are commonly referred to as hybrid meta-heuristics [12]. The surveys on hybrid algorithms that combine the PSO and differential evolution (DE) [13] were presented recently [14, 15]. These PSO-DE hybrids usually employ DE to adjust the particle position. But the convergence performance is dependent on the particle velocity. Limiting the velocity can help the particle to get out of local optima traps [16, 17]. In this paper, we will combine these two optimization algorithms and propose the novel hybrid algorithm H-PSO-DE. The DE is employed to regulate the particle velocity rather than the traditional particle position in case that the optimal result has not improved after several iterations. The hybrid algorithm aims to aggregate the advantages of both algorithms to efficiently tackle the optimization problem.

The remainder of this paper is organized as follows. Section 2 briefly describes the basic operations of the PSO, H-PSO, and DE algorithms. Section 3 presents an analysis of the relationship of particle velocity and convergence. Section 4 provides the hybrid optimization method: H-PSO-DE. Section 5 reveals the simulations and analysis of H-PSO-DE in solving unconstrained optimization problems. Finally, conclusions are given in Section 6.

2. The PSO, H-PSO, and DE Algorithms

2.1. The PSO Algorithm

The PSO [18–20] is a stochastic population-based optimization approach. Each particle is a 𝐷-dimensional vector, and it consists of a position vector 𝑥_𝑛, which represents a candidate solution of the optimization problem, a velocity vector 𝑣_𝑛, and a memory vector 𝑦_𝑛, which is the best candidate solution encountered by the particle. The velocity and position of the particle are updated in every dimension 𝑑(1⩽𝑑⩽𝐷) by 𝑣_𝑛,𝑑(𝑡+1)=𝑤𝑣_𝑛,𝑑(𝑡)+𝑐₁𝑟₁𝑦_𝑛,𝑑(𝑡)−𝑥_𝑛,𝑑(𝑡)+𝑐₂𝑟₂𝑦_{𝑝(𝑛),𝑑}(𝑡)−𝑥_𝑛,𝑑(𝑡),𝑥_𝑛,𝑑(𝑡+1)=𝑥_𝑛,𝑑(𝑡)+𝑣_𝑛,𝑑(𝑡+1), (1)(2) where 𝑤 is the inertia weight, which determines how much of the previous velocity the particle is preserved. 𝑐₁ and 𝑐₂ are positive constants. 𝑟₁ and 𝑟₂ are randomly chosen numbers uniformly distributed in the interval [0,1]. 𝑦_𝑝(𝑛) represents the best position achieved by any member of the population.

2.2. The H-PSO Algorithm

In H-PSO [21], all particles are arranged in a hierarchy. The hierarchy is defined by the height h, the branching degree bd, and the total number of nodes tnn of the corresponding tree.

In H-PSO, the iteration starts with the evaluation of the objective function of each particle at its current position. Then, the new velocity vectors and the new positions for the particles are determined. This means that for particle 𝑛, the value of 𝑦_𝑝(𝑛) in (1) equals 𝑦_𝑚, with 𝑚 being the particle in the parent node of the node of particle 𝑛. H-PSO uses 𝑦_𝑝(𝑛)=𝑦_𝑛 only when particle 𝑛 is in the root. If the function value of a particle 𝑛 is better than the function value at its personal best position so far, then the new position is stored in 𝑦_𝑛. For each particle 𝑛 in a node of the tree, its own best solution is compared to the best solution found by the particles in the child nodes 𝑆(𝑛). If the best of these particles 𝑚 is better than particle 𝑛, then particles 𝑛 and 𝑚 swap their places within the hierarchy.

2.3. The DE Algorithm

The DE [11, 13, 22] is a stochastic parallel direct search method. More specifically, DE’s basic strategy can be summarized as follows.

Initialization. DE begins with a randomly initiated population of 𝑁  𝐷-dimensional parameter vectors 𝑥_𝑖,𝑔, 𝑖=1,2,…,𝑁 as a population for each generation 𝑔. The initial population (𝑔=0) of the 𝑗th parameter of the 𝑖th vector is 𝑥_𝑗,𝑖,0=𝑥_𝑗,min+rand_𝑖,𝑗[0,1]⋅𝑥_𝑗,max−𝑥_𝑗,min, (3) where 𝑥_𝑗,min and 𝑥_𝑗,max indicate the lower and upper bounds, respectively. rand_𝑖,𝑗[0,1] is a uniformly distributed random number lying between 0 and 1.

Mutation. DE mutates and recombines the population to produce a population of 𝑁 trial vectors. Specifically, for each individual 𝑥_𝑖,𝑔, a mutant vector 𝜐_𝑖,𝑔 is generated according to 𝜐_𝑖,𝑔=𝑥_𝑟^𝑖_1,𝑔+𝐹⋅𝑥_𝑟^𝑖_2,𝑔−𝑥_𝑟^𝑖_3,𝑔, (4) where 𝐹, commonly known as scale factor, is a positive real number. Three other random individuals 𝑥_𝑟^𝑖_1,𝑔, 𝑥_𝑟^𝑖_2,𝑔, and 𝑥_𝑟^𝑖_3,𝑔 are sampled randomly from the current population such that 𝑟^𝑖₁,𝑟^𝑖₂,𝑟^𝑖₃∈{1,2,…,𝑁}, and 𝑖≠𝑟^𝑖₁≠𝑟^𝑖₂≠𝑟^𝑖₃.

Crossover. DE crosses each vector with a mutant vector: 𝑢_{𝑗,𝑖,𝑔}=𝜐_{𝑗,𝑖,𝑔},ifrand_𝑖,𝑗[0,1]≤𝐶_𝑟or𝑗=𝑗_rand,𝑥_{𝑗,𝑖,𝑔},otherwise, (5) where 𝐶_𝑟 is called the crossover rate.

Selection. To decide whether or not it should become a member of generation 𝑔+1, the trial vector 𝜐_𝑖,𝑔 is compared to the target vector 𝑥_𝑖,𝑔 using the greedy criterion. The selection operation is described as 𝑥_𝑖,𝑔+1=𝑢_𝑖,𝑔,if𝑓𝑢_𝑖,𝑔≤𝑓𝑥_𝑖,𝑔,𝑥_𝑖,𝑔,otherwise, (6) where 𝑓(𝑥) is the objective function to be minimized.

3. Relationship of Particle Velocity and Convergence

This section presents an analysis of the relationship of particle velocity and convergence.

Substituting (1) into (2) results in 𝑥_𝑛,𝑑(𝑡+1)=𝑥_𝑛,𝑑(𝑡)+𝑤𝑣_𝑛,𝑑(𝑡)+𝑐₁𝑟₁𝑦_𝑛,𝑑(𝑡)−𝑥_𝑛,𝑑(𝑡)+𝑐₂𝑟₂𝑦_{𝑝(𝑛),𝑑}(𝑡)−𝑥_𝑛,𝑑(𝑡). (7) From (2), it is known that 𝑣_𝑛,𝑑(𝑡)=𝑥_𝑛,𝑑(𝑡)−𝑥_𝑛,𝑑(𝑡−1). (8) Substituting (8) into (7) results in 𝑥_𝑛,𝑑(𝑡+1)=1+𝑤−𝑐₁𝑟₁−𝑐₂𝑟₂𝑥_𝑛,𝑑(𝑡)−𝑤𝑥_𝑛,𝑑(𝑡−1)+𝑐₁𝑟₁𝑦_𝑛,𝑑(𝑡)+𝑐₂𝑟₂𝑦_{𝑝(𝑛),𝑑}(𝑡). (9)

This recurrence relation can be written as a matrix-vector product, so that 𝑥_𝑛,𝑑()𝑥𝑡+1_𝑛,𝑑(𝑡)1=1+𝑤−𝑐₁𝑟₁−𝑐₂𝑟₂−𝑤𝑐₁𝑟₁𝑦_𝑛,𝑑(𝑡)+𝑐₂𝑟₂𝑦_{𝑝(𝑛),𝑑}(𝑡)×𝑥100001_𝑛,𝑑(𝑡)𝑥_𝑛,𝑑()1𝑡−1. (10)

The characteristic polynomial of the matrix in (10) is (1−𝜆)(𝑤−𝜆(1+𝑤−𝑐₁𝑟₁−𝑐₂𝑟₂)+𝜆²), which has a trivial root of 𝜆=1 and two other solutions 𝛼=1+𝑤−𝑐₁𝑟₁−𝑐₂𝑟₂+𝛾2,𝛽=1+𝑤−𝑐₁𝑟₁−𝑐₂𝑟₂−𝛾2, (11) where 𝛾=(1+𝑤−𝑐₁𝑟₁−𝑐₂𝑟₂)²−4𝑤.

Note that 𝛼 and 𝛽 are both eigenvalues of the matrix in (10). The explicit form of the recurrence relation (9) is then given by 𝑥_𝑛,𝑑(𝑡)=𝑘₁+𝑘₂𝛼^𝑡+𝑘₃𝛽^𝑡, (12) where 𝑘₁, 𝑘₂, and 𝑘₃ are constants determined by the initial conditions of the system.

Substituting (12) into (8) results in 𝑣_𝑛,𝑑(𝑡)=ℎ₁𝛼^𝑡+ℎ₂𝛽^𝑡, (13)

where ℎ₁=𝑘₂(1−1/𝛼), ℎ₂=𝑘₃(1−1/𝛽).

Consider lim_𝑡→∞𝑣_𝑛,𝑑(𝑡)=lim_𝑡→∞ℎ₁𝛼^𝑡+ℎ₂𝛽^𝑡,lim_𝑡→∞𝑣_𝑛,𝑑(𝑡)=‖𝛼‖,𝛽ℎ0,ifmax<1,₁orℎ₂orℎ₁+ℎ₂,ifmax‖𝛼‖,𝛽=1. (14)(15)

Equation (15) implies that if the PSO algorithm is convergent, the velocity of the particles will decrease to zero or stay unchanged until the end of the iteration.

4. The Proposed H-PSO-DE Algorithm

The main idea of the hybrid H-PSO-DE algorithm is to employ the DE to regulate the particle velocity rather than the traditional particle position in case that the optimal result has not improved after several iterations. If the swarm is going to be in equilibrium, the evolution process will be stagnated as time goes on. To prevent the trend, if the stagnating step of evolution process 𝑔₀ is larger than threshold value 𝐺₀, the particle velocity performs mutation operators. The velocity and position of the particles are updated as follows.

If (rand()<𝐶_𝑟 or 𝑑==𝑘, 𝑘∈[1,𝐷]), then 𝑣^′_𝑛,𝑑=𝑣_𝑛,𝑑+𝐹⋅𝑣_1,𝑑−𝑣_2,𝑑,𝑥^′_𝑛,𝑑=𝑥_𝑛,𝑑+𝑣^′_𝑛,𝑑, (16) where 𝑥_𝑛,𝑑=𝑟⋅𝑦_𝑛,𝑑+(1−𝑟)⋅𝑦_{𝑝(𝑛),𝑑}, 𝑟 is a random number in the interval [0,1], and 𝑣_1,𝑑 and 𝑣_2,𝑑 are sampled randomly from 𝑣_𝑛.

The procedure for H-PSO-DE algorithm is presented in Algorithm 1.

Step 1. Initialize particles in swarm. Initialize location 𝑥𝑛, and velocity 𝑣𝑛    of each particle 𝑛. Best position 𝑦𝑛=𝑥𝑛. Step 2. Evaluate objective function 𝑓(𝑥𝑛) and update personal best.    for each particle 𝑛, do       if 𝑓(𝑥𝑛)<𝑓(𝑦𝑛) then      𝑦𝑛=𝑥𝑛     end if    end for Step 3. Swap particles.      for each particle 𝑛, do      Determine the best successor 𝑚= arg      min 𝑓𝑦𝑞∣𝑞∈𝑆(𝑛), where are the successors of 𝑛.     if 𝑓(𝑦𝑚)<𝑓(𝑦𝑛) then      Swap particles 𝑛 and 𝑚     end if    end for Step 4. Update the position and velocity of the 𝑛th particle.    for each particle 𝑛, do   Update the velocity 𝑣𝑛 in each dimension 𝑑:     𝑣𝑛,𝑑=𝑤𝑣𝑛,𝑑+𝑐1𝑟1𝑦𝑛,𝑑−𝑥𝑛,𝑑+𝑐2𝑟2𝑦𝑝(𝑛),𝑑−𝑥𝑛,𝑑     Move the particle: 𝑥𝑛,𝑑=𝑥𝑛,𝑑+𝑣𝑛,𝑑    end for Step 5. Judge the evolution process of H-PSO.    if (𝑔0==𝐺0), then     goto Step 6     else goto Step 7    end if Step 6. Update the velocity and position of the particle according to (16).    if 𝑓(𝑥′𝑛)<𝑓(𝑥𝑛) then 𝑥𝑛=𝑥′𝑛    end if Step 7. If a stopping criterion is met, then output the global best position and stop;    otherwise, repeat Step 2–Step 6.

Algorithm 1 

Procedure for the H-PSO-DE.

5. Simulations and Results

In this section, we present a simulation study to validate the proposed H-PSO-DE algorithm. A set of test functions that are commonly used in the field of continuous function optimization is listed in the appendix. They are a set of curvilinear functions for difficult unconstrained minimization problems. For illustration, the landscapes of two-dimensional versions of the six functions are depicted in Figure 1. The first two functions (Sphere and Rosenbrock) are unimodal functions, and they have a single local optimum that is also the global optimum. The remaining functions are multimodal, and they have several local optima. Note that the dimensional increase of these scalable functions does not change their basic features.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 1 

An illustration for 2-dimensional landscapes of the test functions. (a) Sphere function; (b) Rosenbrock function; (c) Rastrigin function; (d) Griewank function; (e) Ackley function; and (f) Schaffer’s F6.

In our experiments, the H-PSO uses the parameter values 𝑤=0.729, and 𝑐₁=𝑐₂=1.494 as suggested in [23] for a faster convergence rate. The population size that has been used is 𝑁=21. The maximal number of generations uses 𝐺=5000. The remainder parameters are set as 𝐶_𝑟=0.5, 𝐹=0.6, 𝑟=0.5, 𝐺₀=8, ℎ=3, and 𝑏𝑑=4. Thirty independent runs were carried out. The convergence behavior of the H-PSO is shown in Figure 2. For comparison purpose, the H-PSO-DE is also given in the same figure. As shown in Figure 2, the convergence performance of the H-PSO-DE is better than the H-PSO. H-PSO-DE is compared with H-PSO, DE, and PSO-DE [1] in terms of the selected performance metrics, such as the mean, maximum, and minimum values. In DE, we use DE/rand/1/bin strategy (𝐶_𝑟=0.5, 𝐹=0.6). As shown in Tables 1, 2, and 3, the H-PSO-DE outperforms H-PSO, DE, and PSO-DE. The H-PSO-DE is quite competitive when compared with the other existing methods.

Table 1 

Comparing the mean value of H-PSO-DE with respect to the other state-of-the-art algorithms.

Table 2 

Comparing the maximum value of H-PSO-DE with respect to the other state-of-the-art algorithms.

Table 3 

Comparing the minimum value of H-PSO-DE with respect to the other state-of-the-art algorithms.

Figure 2 

Convergence graph of the H-PSO and the H-PSO-DE for 𝑓1–𝑓6.

6. Conclusions

In this paper, a new method named H-PSO-DE is proposed to solve optimization problems, which improves the performance of the H-PSO by incorporating DE. In H-PSO-DE, when the evolution process is stagnated for several generations, all the particles may lose the ability of finding a better solution. Then, the DE is employed to regulate the particle velocity to avoid wasting too much calculation time for vain search, so the searching efficiency of the H-PSO-DE is improved greatly. The H-PSO-DE is compared on test functions with H-PSO, DE, and PSO-DE. It is shown that H-PSO-DE performs significantly better.

Appendix

Benchmark Functions

Sphere: 𝑓₁(𝑥)=³⁰_𝑖=1𝑥²_𝑖,−50≤𝑥_𝑖𝑓≤50,min₁=𝑓₁(0,…,0)=0. (A.1)

Rosenbrock: 𝑓₂(𝑥)=²⁹_𝑖=1100𝑥_𝑖+1−𝑥²_𝑖²+𝑥_𝑖−1²,−30≤𝑥_𝑖𝑓≤30,min₂=𝑓₂(1,…,1)=0. (A.2)

Rastrigin: 𝑓₃(𝑥)=³⁰_𝑖=1𝑥²_𝑖−10cos2𝜋𝑥_𝑖+10,−5.12≤𝑥_𝑖𝑓≤5.12,min₃=𝑓₃(0,…,0)=0. (A.3)

Griewank: 𝑓₄(𝑥)=14000³⁰_𝑖=1𝑥²_𝑖−Π³⁰_𝑖=1𝑥cos^𝑖√𝑖+1,−600≤𝑥_𝑖𝑓≤600,min₄=𝑓₄(0,…,0)=0. (A.4)

Ackley: 𝑓₅(𝑥)=−20exp−0.2130³⁰_𝑖=1𝑥²_𝑖1−exp30³⁰_𝑖=1cos2𝜋𝑥_𝑖−32≤𝑥_𝑖𝑓≤32,min₅=𝑓₅(0,…,0)=0. (A.5)

Schaffer’s F6: 𝑓₆(𝑥)=0.5+sin𝑥²₁+𝑥²₂−0.51+0.001𝑥²₁+𝑥²₂²,−100≤𝑥_𝑖𝑓≤100,min₆=𝑓₆(0,0)=0. (A.6)

Acknowledgments

This work was supported by the Key Project of Chinese Ministry of Education (no. 212135), the Guangxi Natural Science Foundation (no. 2012GXNSFBA053165), the project of Education Department of Guangxi (no. 201203YB131), and the Doctoral Initiating Project of Guangxi University of Science and Technology (no. 11Z09).