# Training ANFIS Model with an Improved Quantum-Behaved Particle Swarm Optimization Algorithm

Authors: Peilin Liu; Wenhao Leng; Wei Fang
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

This paper proposes a novel method of training the parameters of adaptive-network-based fuzzy inference system (ANFIS). Different from the previous works which emphasized on gradient descent (GD) method, we present an approach to train the parameters of ANFIS by using an improved version of quantum-behaved particle swarm optimization (QPSO). This novel variant of QPSO employs an adaptive dynamical controlling method for the contraction-expansion (CE) coefficient which is the most influential algorithmic parameter for the performance of the QPSO algorithm. The ANFIS trained by the proposed QPSO with adaptive dynamical CE coefficient (QPSO-ADCEC) is applied to five example systems. The simulation results show that the ANFIS-QPSO-ADCEC method performs much better than the original ANFIS, ANFIS-PSO, and ANFIS-QPSO methods.


---

## Body

## 1. Introduction

Fuzzy systems (FSs) have been successfully applied in many areas, such as system modeling and controls. To ease the design and improve system performance, many neural or statistical learning approaches that automatically generate fuzzy rules have been proposed [1]. A fuzzy inference system employing fuzzy if-then rules can model the qualitative aspects of human knowledge and reasoning processes without using precise quantitative analyses. This fuzzy modeling or fuzzy identification, first explored systematically by Sugeno and Kang [2], has found numerous practical applications in control [3, 4], prediction, and inference [5, 6]. However, there are some basic aspects of this approach which are in need of better understanding. First, there are no standard methods for transforming human knowledge or experience into the rule base and data base of a fuzzy inference system. Second, there is a need for effective methods for tuning the membership functions (MF) so as to minimize the output error measure or maximize performance index.This paper concentrates on a so-called adaptive-network-based fuzzy inference system (ANFIS), which can serve as a basis for constructing a set of fuzzy if-then rules with appropriate membership functions to generate the stipulated input-output pairs [7]. The ANFIS model is a reprehensive of adaptive fuzzy systems which are generated through training processes and are known as evolving fuzzy systems. It has been one of attractive researches focusing on fuzzy systems in recent years.The TSK [2] is a fuzzy system with crisp functions and has been found to be efficient in complex applications [4]. It has been proved that, with proper number of rules, a TSK system is able to approximate every plant. As such, TSK systems are widely used in the ANFIS and play the advantage of good applicability since they can be interpreted as local linearization modeling and conventional linear techniques for state estimation and control.The ANFIS has both the advantages of neural networks and fuzzy systems. However, training the parameters of the ANFIS model is one of the main issues encountered when the model is applied to the real-world problems. Most of the training methods for the ANFIS are based on gradient descent (GD) approaches, where calculation of gradient in each step is tractable since the chain rule used may cause many local minima of the problem. The gradient methods are known to be local search approaches and their performances generally depend on initial values of parameters so that it is difficult for them to find the global optimal model parameters. Since the design of FSs can be reduced to an optimization problem, many researchers have proposed to design FSs by employing metaheuristics such as genetic algorithms (GAs) [8–10] and particle swarm optimization (PSO) [11–13]. However, GAs have always been complaint about their slow convergence speed, while PSO may encounter premature convergence at the later stage of the search process and is sensitive to neighborhood topology.This paper explores the applicability of a variant of PSO, quantum-behaved particle swarm optimization (QPSO), to training of the ANFIS model. The QPSO algorithm was inspired by quantum mechanics [14–17]. Its iterative equation is very different from that of PSO and can lead QPSO to be globally convergent [18–20]. Besides, unlike PSO, QPSO needs no velocity vectors for particles and has fewer parameters to adjust, making it easier to implement. The QPSO algorithm has been shown to successfully solve a wide range of continuous optimization problems. Many empirical studies show that QPSO has stronger global search ability when solving various continuous optimization problems [21–27]. In order to make a further improvement of the QPSO, in this paper, we propose an adaptive dynamical control method for the contraction-expansion (CE) coefficient of the algorithm, which is the most influential algorithmic parameter that can be tuned to adjust the balance between the local search and global search of the particle. The improved QPSO algorithm is applied to ANFIS training and is tested on several example systems, with performance comparison between the improved QPSO, the original QPSO, and PSO algorithms.The rest of the paper is organized as follows: in Section2, we present a review of fuzzy systems and the ANFIS model. Section 3 presents the principle of QPSO. The improved QPSO approach is proposed in Section 4. Application of the QPSO and its improved version to ANFIS model training is described in Section 5. Section 6 presents the simulation results for five example systems using ANFIS model trained by the optimization algorithms. Finally, the paper is concluded in Section 7.


## 2. Fuzzy Systems and ANFIS Model

### 2.1. Fuzzy Systems

This subsection describes the fuzzy systems to be optimized through training processes in this study and their mathematical expressions. The TSK fuzzy model was originally proposed by Sugeno and Kang in an effort to formalize a systematic approach to generating fuzzy rules from an input-output data set [2]. The Takagi-Sugeno-Kang (TSK) fuzzy system can be of zero order or first order. The ith rule, which is denoted Ri, in a TSK fuzzy system is represented in the following form: Ri
: ifx1(k) is Ai1 and⋯and xn(k) is Ain, then y(k) is fi, i=1,…,r,where k is the number of time step, x1(k),…,xn(k) are input variables, y(k) is the output variable of the system, Aij is a fuzzy set, and fi is the consequent part function. The fuzzy set Aij uses a bell-shaped membership function given by
(1)μ(x)=11+((xj-cij)/aij)2bij,
where ai, bi, and ci are the parameter set of the fuzzy set Aij. Each of these parameters has a physical meaning: ci determines the center of the corresponding membership function,  ai is the half width, and bi (together with ai) controls the slopes at the crossover points. In the inference engine, the fuzzy AND operation is implemented by the algebraic product in fuzzy theory. Thus, given an input dataset x=(x1,…,xn), the firing strength μ(x) of rule i is calculated by
(2)μ(x¯)=∏j=1nMij(x)=∏j=1n{11+((xj-cij)/aij)2bij}.For zero-order and first-order TSK-type fuzzy systems, the consequent functionfi is set to a real value pi0 and a linear function of input variables pi0+∑j=1npijxj for each i, respectively. If there are r rules in a fuzzy system, the output of the system, which is calculated by the weighted-average-defuzzification method, is given by
(3)y=∑i=1rμi(x¯)fi∑i=1rμi(x)i.In this work, the proposed improved QPSO algorithm and other optimization algorithms will be used to optimize free parametersai,  bi, and ci in each rule.


### 2.2. The Architecture of ANFIS

This subsection presents the architecture of ANFIS network. A detailed coverage of ANFIS can be found in [7]. The ANFIS network is a neurofuzzy network that was proposed by Jang in 1993 [6]. Since the ANFIS is an adaptive network, parts of its nodes are adaptive, which means that their outputs depend on the parameters belonging to these nodes. Two kinds of learning algorithms have been proposed to tune these parameters to optimize the approaching performance during the training period. For simplicity, the above-mentioned system is supposed to have two inputs and one output, and its rule base contains two fuzzy if-then rules of TSK fuzzy model. A typical TSK fuzzy model has the two rules that can be stated asRule 1: ifx is A1 and y is B1, then f1=p1x+q1y+r1,Rule 2: ifx is A2 and y is B2, then f2=p2x+q2y+r2,where x and y are the inputs of the ANFIS, A and B are the fuzzy sets, and fi (i=1,2) is a first-order polynomial and represents the outputs of the first-order TSK fuzzy inference system. In the above rules, pi, qi, and ri (i=1,2) are the parameters set, referred to as the consequent parameters.The architecture of ANFIS is shown in Figure1, and the node function in each layer is described below.Figure 1
The architecture of ANFIS network.Layer 1. This layer is the layer of membership functions that contains adaptive nodes with node functions described as
(4)Oi1=μAi(x),Oi1=μB(i-2)(y),(i=1,2),
where x and y are the input nodes, A and B are the linguistic labels, μ(x) and μ(y) are the membership functions which usually adopts a bell shape with the maximum and minimum values equal to 1 and 0, respectively:
(5)μ(x)=11+((x-ci)/ai)2bi,
where ai, bi, and ci are the parameters set. When values of these parameters change, the bell-shaped functions vary accordingly, exhibiting various forms of membership functions on linguistic labels A and B. In fact, any continuous and piecewise differentiable functions, such as commonly used trapezoidal or triangular-shaped membership functions, are also qualified candidates for node functions in this layer. Parameters in this layer are referred to as premise parameters.Layer 2. Every node in this layer is a fixed node, marked by a circle and labeled Π, with the node function to be multiplied by input signals to serve as output. Consider
(6)Oi2=ωi=μAi(x)·μBi(y),(i=1,2).
The output ωi represents the firing strength of a rule. The output of each node represents the firing strength of a rule. In fact, other T-norm operators that perform the generalized AND can be used as the node function in this layer.Layer 3. Every node in this layer is a fixed node, marked by a circle and labeled N, with the node function to normalize the firing strength by calculating the ratio of the ith node firing strength to the sum of all rules’ firing strength. Moreover,
(7)Oi3=ϖi=ωi∑ωi=ω1ω1+ω2,(i=1,2).
For convenience, the outputs of this layer are referred to as normalized firing strengths.Layer 4.Every node in this layer is an adaptive node, marked by a square, with node function given by
(8)Oi4=ϖi·fi=ϖi(pix+qiy+ri),(i=1,2),
where ϖ  is the output of  layer 3 and  {pi,qi,ri} is the parameter set. The parameters in this layer are referred to as consequent parameters.Layer 5.Every node in this layer is a fixed node, and the overall output can be expressed as linear combination of the consequent parameters. Consider
(9)Oi5=fout=∑iϖi·fi=ϖ1f1+ϖ2f2=ω1ω1+ω2f1=ω2ω1+ω2f2=(ϖ1x)p1+(ϖ1y)q1+(ϖ1)r1+(ϖ2x)p2+(ϖ2y)q2+(ϖ2)r2.


### 2.3. Hybrid Learning Algorithm for ANFIS Model

It can be seen that there are two modifiable parameter sets,{ai,bi,ci} labeled as premise parameters and  {pi,qi,ri} labeled as consequent parameters. The aim of the training procedure for this architecture is to tune the above two parameter sets to make the ANFIS output fit the training data. Each epoch of this hybrid learning procedure is composed of two passes: a forward pass and a backward pass. In the forward pass, the premise parameters are fixed, and the least squares estimation (LSE) is applied to identify consequent parameters. When the optimal parameters are found, the backward pass starts with the consequent parameters fixed, the error rate of output node back propagates from output end toward the input end, and the premise parameters are updated by the gradient descent (GD) method.These methods update premise parameters by using GD or Kalman filtering and appear to be prone to trap into the local optima. In this paper, we propose a QPSO with adaptive dynamical contraction-expansion coefficient (QPSO-ADCEC) and employ this algorithm to train the parameters of ANFIS for the purpose of obtaining the global optimal solution.


## 2.1. Fuzzy Systems

This subsection describes the fuzzy systems to be optimized through training processes in this study and their mathematical expressions. The TSK fuzzy model was originally proposed by Sugeno and Kang in an effort to formalize a systematic approach to generating fuzzy rules from an input-output data set [2]. The Takagi-Sugeno-Kang (TSK) fuzzy system can be of zero order or first order. The ith rule, which is denoted Ri, in a TSK fuzzy system is represented in the following form: Ri
: ifx1(k) is Ai1 and⋯and xn(k) is Ain, then y(k) is fi, i=1,…,r,where k is the number of time step, x1(k),…,xn(k) are input variables, y(k) is the output variable of the system, Aij is a fuzzy set, and fi is the consequent part function. The fuzzy set Aij uses a bell-shaped membership function given by
(1)μ(x)=11+((xj-cij)/aij)2bij,
where ai, bi, and ci are the parameter set of the fuzzy set Aij. Each of these parameters has a physical meaning: ci determines the center of the corresponding membership function,  ai is the half width, and bi (together with ai) controls the slopes at the crossover points. In the inference engine, the fuzzy AND operation is implemented by the algebraic product in fuzzy theory. Thus, given an input dataset x=(x1,…,xn), the firing strength μ(x) of rule i is calculated by
(2)μ(x¯)=∏j=1nMij(x)=∏j=1n{11+((xj-cij)/aij)2bij}.For zero-order and first-order TSK-type fuzzy systems, the consequent functionfi is set to a real value pi0 and a linear function of input variables pi0+∑j=1npijxj for each i, respectively. If there are r rules in a fuzzy system, the output of the system, which is calculated by the weighted-average-defuzzification method, is given by
(3)y=∑i=1rμi(x¯)fi∑i=1rμi(x)i.In this work, the proposed improved QPSO algorithm and other optimization algorithms will be used to optimize free parametersai,  bi, and ci in each rule.


## 2.2. The Architecture of ANFIS

This subsection presents the architecture of ANFIS network. A detailed coverage of ANFIS can be found in [7]. The ANFIS network is a neurofuzzy network that was proposed by Jang in 1993 [6]. Since the ANFIS is an adaptive network, parts of its nodes are adaptive, which means that their outputs depend on the parameters belonging to these nodes. Two kinds of learning algorithms have been proposed to tune these parameters to optimize the approaching performance during the training period. For simplicity, the above-mentioned system is supposed to have two inputs and one output, and its rule base contains two fuzzy if-then rules of TSK fuzzy model. A typical TSK fuzzy model has the two rules that can be stated asRule 1: ifx is A1 and y is B1, then f1=p1x+q1y+r1,Rule 2: ifx is A2 and y is B2, then f2=p2x+q2y+r2,where x and y are the inputs of the ANFIS, A and B are the fuzzy sets, and fi (i=1,2) is a first-order polynomial and represents the outputs of the first-order TSK fuzzy inference system. In the above rules, pi, qi, and ri (i=1,2) are the parameters set, referred to as the consequent parameters.The architecture of ANFIS is shown in Figure1, and the node function in each layer is described below.Figure 1
The architecture of ANFIS network.Layer 1. This layer is the layer of membership functions that contains adaptive nodes with node functions described as
(4)Oi1=μAi(x),Oi1=μB(i-2)(y),(i=1,2),
where x and y are the input nodes, A and B are the linguistic labels, μ(x) and μ(y) are the membership functions which usually adopts a bell shape with the maximum and minimum values equal to 1 and 0, respectively:
(5)μ(x)=11+((x-ci)/ai)2bi,
where ai, bi, and ci are the parameters set. When values of these parameters change, the bell-shaped functions vary accordingly, exhibiting various forms of membership functions on linguistic labels A and B. In fact, any continuous and piecewise differentiable functions, such as commonly used trapezoidal or triangular-shaped membership functions, are also qualified candidates for node functions in this layer. Parameters in this layer are referred to as premise parameters.Layer 2. Every node in this layer is a fixed node, marked by a circle and labeled Π, with the node function to be multiplied by input signals to serve as output. Consider
(6)Oi2=ωi=μAi(x)·μBi(y),(i=1,2).
The output ωi represents the firing strength of a rule. The output of each node represents the firing strength of a rule. In fact, other T-norm operators that perform the generalized AND can be used as the node function in this layer.Layer 3. Every node in this layer is a fixed node, marked by a circle and labeled N, with the node function to normalize the firing strength by calculating the ratio of the ith node firing strength to the sum of all rules’ firing strength. Moreover,
(7)Oi3=ϖi=ωi∑ωi=ω1ω1+ω2,(i=1,2).
For convenience, the outputs of this layer are referred to as normalized firing strengths.Layer 4.Every node in this layer is an adaptive node, marked by a square, with node function given by
(8)Oi4=ϖi·fi=ϖi(pix+qiy+ri),(i=1,2),
where ϖ  is the output of  layer 3 and  {pi,qi,ri} is the parameter set. The parameters in this layer are referred to as consequent parameters.Layer 5.Every node in this layer is a fixed node, and the overall output can be expressed as linear combination of the consequent parameters. Consider
(9)Oi5=fout=∑iϖi·fi=ϖ1f1+ϖ2f2=ω1ω1+ω2f1=ω2ω1+ω2f2=(ϖ1x)p1+(ϖ1y)q1+(ϖ1)r1+(ϖ2x)p2+(ϖ2y)q2+(ϖ2)r2.


## 2.3. Hybrid Learning Algorithm for ANFIS Model

It can be seen that there are two modifiable parameter sets,{ai,bi,ci} labeled as premise parameters and  {pi,qi,ri} labeled as consequent parameters. The aim of the training procedure for this architecture is to tune the above two parameter sets to make the ANFIS output fit the training data. Each epoch of this hybrid learning procedure is composed of two passes: a forward pass and a backward pass. In the forward pass, the premise parameters are fixed, and the least squares estimation (LSE) is applied to identify consequent parameters. When the optimal parameters are found, the backward pass starts with the consequent parameters fixed, the error rate of output node back propagates from output end toward the input end, and the premise parameters are updated by the gradient descent (GD) method.These methods update premise parameters by using GD or Kalman filtering and appear to be prone to trap into the local optima. In this paper, we propose a QPSO with adaptive dynamical contraction-expansion coefficient (QPSO-ADCEC) and employ this algorithm to train the parameters of ANFIS for the purpose of obtaining the global optimal solution.


## 3. Quantum-Behaved Particle Swarm Optimization

In the PSO withM individuals, each individual is treated as a volume-less particle in the D-dimensional space, with the current position vector and velocity vector of particle i at the nth iteration represented as Xi,n=(Xi,n1,Xi,n2,…,Xi,nD) and Vi,n=(Vi,n1,Vi,n2,…,Vi,nD) [28, 29]. The particle moves according to the following equations:
(10)Vi,n+1j=wVi,nj+c1ri,nj(Pi,nj-Xi,nj)+c2Ri,nj(Gi,nj-Xi,nj),(11)Xi,n+1j=Xi,nj+Vi,n+1j,
for i=1,2,.…,M and j=1,2,…,D, where c1 and c2 are known as acceleration coefficients. The parameter w is known as the inertia weight which can be adjusted to balance the explorative search and the exploitive search of the particle. Vector Pi,n=(Pi,n1,Pi,n2⋯Pi,nD) is the best previous position (the position giving the best objective function value or fitness value) of particle i and called personal best (pbest) position, and vector Gn=(Gn1,Gn2⋯GnD) is the position of the best particle among all the particles in the population and called global best (gbest) position. Without loss of generality, we consider the following maximization problem:
(12)Minimizef(X)s.t.X∈S⊆RD,
where f(X) is an objective function continuous almost everywhere and S is the feasible space. Accordingly, Pi,n can be updated by
(13)Pi,n={Xi,niff(Xi,n)<f(Pi,n-1),Pi,n-1iff(Xi,n)≥f(Pi,n-1),
and Gn can be found by Gn=Pg,n, where g=argmax1≤i≤M[f(pi,n)]. The parameters ri,nj and  Ri,nj are sequences of two different sequences of random numbers distributed uniformly within (0,1), which is denoted by ri,nj,Ri,nj~U(0,1). Generally, the value of Vi,nj is restricted in the interval [-Vmax,Vmax].Trajectory analysis in [30] showed that convergence of the PSO algorithm may be achieved if each particle converges to its local attractor, pi,n=(pi,n1,pi,n2⋯pi,nD), defined at the coordinates
(14)pi,nj=c1ri,njpi,nj+c2Ri,njGnjc1ri,njc2Ri,nj1≤j≤D
or
(15)pi,nj=φi,nj·pi,nj+(1-φi,nj)·Gnj,
where  φi,nj=c1ri,nj/(c1ri,nj+c2Ri,nj) with regard to the random numbers ri,nj and Ri,nj in (10). In PSO, the acceleration coefficients c1 and c2 are generally set to be equal; that is, c1=c2, and thus φi,nj is a sequence of uniformly distributed random numbers over (0,1). As a result, (15) can be restated as
(16)pi,nj=φi,nj·Pi,nj+(1-φi,nj)·Gnj,ϕi,nj~U(0,1).In QPSO, each single particle is treated as a spin-less one moving in quantum space. Thus state of the particle is characterized by wave functionψ, where |ψ|2 is the probability density function of its position. Inspired by convergence analysis of the particle in PSO [30], it is assume that, at the nth iteration, particle i flies in the D-dimensional space with δ potential well centered at pi,nj on the jth dimension  (1≤j≤D). Let  Yi,n+1j=|Xi,nj-pi,nj|; we can obtain the normalized wave function at iteration n+1(17)ψ(Yi,n+1j)=1Li,njexp(-Yi,n+1jLi,nj),
which satisfies the bound condition that ψ(Yi,n+1j)→0 as Yi,n+1j→∞. Li,nj is the characteristic length of the wave function. By the statistical interpretation of wave function, the probability density function is given by
(18)Q(Yi,n+1j)=|ψ(Yi,nj)|2=1Li,njexp(-2Yi,n+1jLi,nj),
and thus the probability distribution function is
(19)F(Yi,n+1j)=1-exp(-2Yi,n+1jLi,nj).
Using Monte Carlo method, we can measure thejth component of position of particle i at the (n+1)th iteration by
(20)Xi,n+1j=pi,n+1j±Li,nj2ln(1ui,n+1j),ui,n+1j~U(0,1),
where ui,n+1j is a sequence of random numbers uniformly distributed over (0,1). The value of Li,nj is determined by
(21)Li,nj=2α·|Xi,nj-Ci,nj|,
where Cn=(Cn1,Cn2,…,CnD) is called mean best (mbest) position defined by the average of the pbest positions of all particles; namely,
(22)Cnj=(1M)∑i=1Mpi,nj(1≤j≤D).
Thus the position of the particle updates according to the following equation:
(23)Xi,n+1j=pi,nj±α·|Xi,nj-Ci,nj|·ln(1ui,n+1j).
The parameter α in (21) and (23) is called contraction-expansion (CE) coefficient, which can be adjusted to balance the local search and the global search of the algorithm during the optimization process. The current position of the particle in QPSO is thus updated according to (16) and (23).The QPSO algorithm starts with the initialization of the particle’s current positions and theirpbest positions (setting Pi,0=Xi,0), followed by the iteration of updating the particle swarm. At each iteration, the mbest position of the particle swarm is computed and the current position of each particle is updated according to (16) and (23). Before each particle updates its current position, its fitness value is evaluated, and then its pbest position and the current gbest position are updated. In (23), the probability of using either operation ‘‘+’’ or operation ‘‘−’’ is equal to 0.5. The search procedure continues until the termination condition is met.We outline the procedure of the QPSO algorithm as follows.


### 3.1. Procedure of the QPSO

Step 1.
Initialize the population; that is, initialize the current position and personal best position of each particle.Step 2.
Execute the following steps.Step 3.
Compute mean best positionCn according to (22).Step 4.
Properly select the value ofα.Step 5.
For each particle in the population, execute from Step6 to Step 8.Step 6.
Evaluate the objective function value of the current position of the particle, that is,f(Xi,n).Step 7.
UpdatePi,n and Gn according to  (13).Step 8.
Update each component of the particle’s position according to (16) and (23).Step 9.
While the termination condition is not met, return to Step2.Step 10.
Output the results.


## 3.1. Procedure of the QPSO

Step 1.
Initialize the population; that is, initialize the current position and personal best position of each particle.Step 2.
Execute the following steps.Step 3.
Compute mean best positionCn according to (22).Step 4.
Properly select the value ofα.Step 5.
For each particle in the population, execute from Step6 to Step 8.Step 6.
Evaluate the objective function value of the current position of the particle, that is,f(Xi,n).Step 7.
UpdatePi,n and Gn according to  (13).Step 8.
Update each component of the particle’s position according to (16) and (23).Step 9.
While the termination condition is not met, return to Step2.Step 10.
Output the results.


## 4. The QPSO with Adaptively Dynamical CE Coefficient (QPSO-ADCEC)

The CE coefficient is the most influential algorithmic parameter for the search performance of the QPSO. In [19], the influence of the CE coefficient on particles’ dynamical behaviour and the algorithmic performance were theoretically and empirically analyzed, and two control methods for the CE coefficient, linearly decreasing and fixed value methods, were investigated in depth on a suite of well-known benchmark functions. Here, we propose a novel control method for the parameter and thus propose a modified QPSO. The motivation of this proposal is to control the CE coefficient in an adaptive and dynamical way according to evolution process and the diversity of the particle swarm during the search process.From experience it is well known that the convergence of the QPSO depends on the evolution speed of the fitness values of the particles. Therefore in order to handle complex and nonlinear optimization process, it would be helpful to include such factor into the CE coefficient design. In this paper, we propose a novel dynamically varying CE coefficient for the QPSO involving two factors, namely, the evolution speed factor and the aggregation degree factor. The evolution speed factor of theith particle at the nth iteration is given by
(24)hi,n=min(f(Pi,n-1),f(Pi,n))max(f(Pi,n-1),f(Pi,n)),
where f(Pi,n) is the fitness value of the personal best position of particle i at iteration n and n is the current iteration number. Note that the evolutionary speed factor lies in (0,1]. This parameter reflects the run time history of the algorithm and the evolution speed of the particle. For the minimization problem described in (12), if the given optimization is a minimization one, f(Pi,n)≤f(Pi,n-1) for each n>0. Thus hi,n can be simply expressed as hi,n=f(Pi,n)/f(Pi,n-1). The smaller the value of hi,n, the faster the decrease of f(Pi,n) and the faster the evolution speed of particle i. After certain number of iterations, the value of hi,n attains its maximum value indicating that f(Pi,n) is not improved any more, and the algorithm stagnates or finds the optimum solution. The aggregation degree factor s of particle i at the nth iteration is defined by
(25)si,n=min(f(Bn),f-n)max(f(Bn),f-n),
where f-n is the average fitness value of all the particles at the nth iteration, that is;
(26)f-n=1M∑i=1Mf(Xi,n),
where M is the swarm size. f(Bn) represents the optimal value that the swarm found in the nth iteration. Note that f(Bn) cannot be replaced by f(Gn) since f(Gn) denotes the optimal value that the entire swarm has found up to the nth iteration. It is obvious that the value of si,n falls in (0,1]. The aggregation can reflect not only the aggregation degree but also the diversity of the swarm in terms of the fitness values of the current position of the particles. The larger the value of si,n, the less the difference between f(Bn) and f-n, and thus the more the aggregation or the smaller the swarm diversity. In the case when all the si,n are equal to 1, all the particles have the same identity with each other in terms of the fitness values of their current positions. The new CE coefficient for each particle can be written as
(27)αi,n=αinitial-(1-hi,n)χ+si,nλ,
where αinitial is the initial value of the CE coefficient and χ and λ typically lie in the range [0,1]. It is obvious that the value of inertia weight lies in between (1-χ) and (1+λ). With this adaptively dynamical CE coefficient, the modified QPSO algorithm, called QPSO with adaptively dynamical CE coefficient (QPSO-ADCEC), can have strong ability to jump out of local search than the original QPSO with linearly decreasing or fixed CE coefficient.


## 5. Learning ANFIS Network by QPSO-ADCEC

This section presents how to employ the QPSO-ADCEC algorithm for learning the ANFIS parameters. The ANFIS has two types of parameters which need to be trained, that is, the premise parameters{ai,bi,ci} and the consequent parameters {pi,qi,ri}. The premise parameters are estimated by the QPSO-ADCEC algorithm, and the consequent parameters are identified by the least squares estimation (LSE). Thus, the position of each particle in the QPSO-ADCEC represents a set of premise parameters. The fitness is defined as root mean squared error (RMSE) between actual output and desired output, which can be expressed by
(28)fitness=∑i=1n(f(i)-f0(i))2n.
The ANFIS coupling with the QPSO-ADCEC algorithm is outlined as below.Step 1.
Initialize the swarm of particles such that the position of each particle are uniformly distributed within the search scope, and set the iteration numbern=0.Step 2.
Set the position of each particleXi,n as the premise parameters {ai,bi,ci}, and identify consequent parameters {pi,qi,ri} with LSE. Then calculate fitness value of each particle, and set each particle’s personal best position as Pi,0=Xi,0.Step 3.
Find out the mean value of all particles’ personal best positionCn  by using (22).Step 4.
For each particle in the population, execute from Step5 to Step 8.Step 5.
Calculate each particle’s fitness value,fitness(Xi,n), and then compare it with the fitness of its personal best position, fitness(Pi,n-1). If fitness(Xi,n)<fitness(Pi,n-1), then Pi,n<Xi,n, otherwise Pi,n<Pi,n-1.Step 6.
Compare the fitness value of each particle’s personal best positionfitness(Pi,n) with that of the global best position, fitness(Gn-1). If fitness(Pi,n)<fitness(Gn-1), then Gn=Pi,n, otherwise Gn=Gn-1.Step 7.
Calculatehi,n and si,n for each particle according to (24) and (25). Then obtain the value of αi,n by using (27).Step 8.
Update the position of each particleXi,n by using (16) and (23), where α is replaced by αi,n.Step 9.
If the termination condition is met, go to exit, otherwise go to Step2, and set n=n+1.


## 6. Simulation Results

To investigate the efficiency of the proposed method, five example systems were tested. The first two examples are the problem of identification of nonlinear systems. For the third example, ANFIS-QPSO-ADCEC was used as an identifier to identify a nonlinear component discrete control system. In the simulation study for the fourth example, we used the proposed method to predict a chaotic time series. The fifth example is to use the proposed ANFIS-QPSO-ADCEC method to predict the carbon dioxide in Box-Jenkins Furnace. For all the examples, besides the ANFIS-QPSO-ADCEC method, the ANFIS model trained with BP algorithm (ANFIS-BP), ANFIS with PSO (ANFIS-PSO), and ANFIS with QPSO (ANFIS-QPSO) were also tested and compared with the proposed ANFIS-QPSO-ADCEC. For the PSO, the inertia weight was set to be 0.73, and the acceleration coefficients were set asc1=c2=1.49. For the QPSO algorithm, the CE coefficient was decreased linearly from 1.0 to 0.5 over the course of search. For the QPSO-ADCEC, the initial CE coefficient αinitial=0.5, the value of χ and λ were set to be 0.6 and 0.4, respectively. These parameters were shown to result in good algorithmic performance in our preliminary experiments. For each of the PSO variants, the swarm size was set to be 50, and the maximum number of iterations was 2000.Example 1 (nonlinear SISO System modeling [1]).
In this example, the nonlinear plant to be identified by first-order TSK-type fuzzy system is described as(29)y=sin(πx)+0.8sin(3πx)+0.2sin(5πx),
where x is the only input. The fuzzy system had five rules with the membership function of each rule being Gaussian function. Each member function had two parameters and there were 10 premise parameters, and 10 consequent parameters in the fuzzy system. Thus the dimensionality of the position of each particle is 10. We chose 100 input data which were randomly generated between −1 and 1, 50 of which were used as the training data and the remainder of which were employed as testing data.
Table1 lists the training error and test error of each method, from which we can see that the training errors (trnRMSE) and test errors (testRMSE) of the ANFIS network trained by the QPSO-ADCEC, QPSO, and PSO are smaller than those of the ANFIS trained by BP algorithm. Among all the PSO variants, the QPSO-ADCEC was able to train the ANFIS model with the smallest training and test errors. Figure 2 visualizes results by using the ANFIS-QPSO-ADCEC for identification of the system and shows that the ANFIS model learned by the QPSO-ADCEC fits the nonlinear system very well.Table 1
Results for Example1.
trnRMSE
testRMSE
ANFIS-BP
1.665
e
-
0
2
1.675
e
-
0
2
ANFIS-PSO
9.796
e
-
0
3
9.905
e
-
0
3
ANFIS-QPSO
6.301
e
-
0
3
7.620
e
-
0
3
ANFIS-QPSO-ADCEC
1.486
e
-
03
2.398
e
-
03Figure 2
The actual output and desired output by ANFIS-QPSO-ADCEC for Example1.Example 2.
Nonlinear MISO System modeling is described as(30)y=sin(πx1)+0.8sin(3πx1)+0.2sin(5πx1)+0.6sin(2πx1)+0.6sin(4πx2)+0.1sin(5πx2)+0.2sin(3πx2)+0.3sin(2πx3)+0.5sin(πx3).
Equation (30) is three-input nonlinear function. Each variable had two membership functions so that the fuzzy system had eight fuzzy rules [1]. Each membership function was Gaussian with 2 parameters, and the fuzzy system had 12 premise parameters and 32 consequent parameters. Thus, the dimensionality of the position of each particle in QPSO-ADCEC was 12. From the grid points of the range [-1,1]×[-1,1]×[-1,1] within the input space of the above equation, 500 data pairs were obtained, where 250 data were used for training and the remaining 250 data for testing.
The results for this example are presented in Table2, where it is shown by training and test errors that the proposed ANFIS method outperformed the other methods on this example system and ANFIS-BP performed worst among all the methods. Figure 3 also shows that the ANFIS model learned by the QPSO-ADCEC fitted the system fairly well.Table 2
Results for Example2.
trnRMSE
testRMSE
ANFIS-BP
4.9757
e
-
03
5.0329
e
-
03
ANFIS-PSO
3.1639
e
-
0
3
3
.
8402
e
-
0
3
ANFIS-QPSO
2
.
9343
e
-
0
3
3
.
239
0
e
-
0
3
ANFIS-QPSO-ADCEC
1.2611
e
-
03
1.8366
e
-
03Figure 3
The actual output and desired output by ANFIS-QPSO-ADCEC for Example2.Example 3 (online identification in control systems [31]).
Here we use ANFIS-QPSO-ADCEC and the other methods to identify a nonlinear component in a control system. The plant under consideration is governed by the following difference equation:(31)y(k+1)=0.3y(k)+0.6y(k-1)+u(k).
We chose y(k) and u(k) as the inputs and y(k+1) as the output. The input to the plant and the model was a sinusoid u(k)=sin(2πk/250), and the adaptation started at k=1 and stopped at k=250. As shown in Figure 4, the output of the model follows the output of the plant almost immediately even after the adaptation stopped at k=250, and the u(k) is changed to u(k)=0.5sin(2πk/250)+0.5sin(2πk/25) after it stopped k=500. Figure 4 also illustrates the results using the ANFIS-QPSO-ADCEC for identification. The training error and test error provided in Table 3 show the superiority of the ANFIS-QPSO-ADCEC to its competitors.Table 3
Results for Example3.
trnRMSE
testRMSE
ANFIS-BP
1
.
22
4
2
e
-
02
1
.
304
1
e
-
02
ANFIS-PSO
5.0413
e
-
03
1.0454
e
-
03
ANFIS-QPSO
3.3281
e
-
0
3
3.3476
e
-
0
3
ANFIS-QPSO-ADCEC
2.1054
e
-
03
2.3770
e
-
03Figure 4
The actual output and desired output by ANFIS-QPSO-ADCEC for Example3.Example 4 (prediction of future values of a chaotic time series [31]).
The simulation results for Examples1 to 3 show that the ANFIS-QPSO-ADCEC can be used to model highly nonlinear functions effectively. In this example, we are to demonstrate how the proposed ANFIS-QPSO-ADCEC can be employed to predict future values of a chaotic time series. The time series used in our simulation is generated by the chaotic Mackey-Glass differential delay equation [22] defined by
(32)x(t)=0.2x(t-τ)1+x10(t-τ)-0.1x(t).
The prediction of future values of this time series is a benchmark problem which has been considered by a number of connectionist researchers (Lapedes and Farber [32], Moody [33], Jones et al. [34], Crowder [35], and Sanger [36]).
The initial conditions forx(0) and τ are 1.2 and 17, respectively. We use four past data for this prediction, and the fuzzy system is generated as
(33)[x(t-18),x(t-12),x(t-6),x(t):x(t+6)].
When t was varying from 118 to 1117, we generated 1000 data pairs for our data and applied the first 500 data pairs for training and the rest 500 data pairs for prediction. Table 4 provides the training error and test error of each method, showing again that the proposed QPSO-ADCEC algorithm performed the best in training the ANFIS model for this example. Figure 5 shows the errors between actual output and desired output by using ANFIS-QPSO-ADCEC.Table 4
Results for Example4.
trnRMSE
testRMSE
ANFIS-BP
2.550
e
-
0
3
2.502
e
-
0
3
ANFIS-PSO
2.073
e
-
0
3
2.043
e
-
0
3
ANFIS-QPSO
1.870
e
-
0
3
1.777
e
-
0
3
ANFIS-QPSO-ADCEC
1.3902
e
-
03
1.3811
e
-
03(a) The actual output and desired output; (b) errors between actual output and desired output by using the ANFIS-QPSO-ADCEC for Example4.
(a)
(b)Example 5 (prediction of carbon dioxide in Box-Jenkins Furnace [33]).
Box-Jenkins Furnace data set contains 296 pairs of input and output data, where the inputx(t) is the flow of methane and the output y(t) is the concentration of carbon dioxide dismissed from the furnace. In our experiment, we chose x(t), x(t-1), x(t-2), y(t-1), and y(t-2) to be the input data and y(t) to be the output of the system [37], thus obtain 294 pairs of input and output data. The first 244 pairs were used to training the ANFIS model and the last 50 pairs were employed to test the trained system.
Table5 presents the training error and test error of each method, showing that the QPSO-ADCEC algorithm performed the best in training the ANFIS model and the trained system had the least test error. Figure 6 plotted the desired output and the actual output of the ANFIS model trained by the QPSO-ADCEC. It can be observed that the obtained system fitted the data very well.Table 5
Results for Example5.
trnRMSE
testRMSE
ANFIS-BP
8.705
e
-
02
1.047
e
+
00
ANFIS-PSO
7.348
e
-
02
1.394
e
-
01
ANFIS-QPSO
3.984
e
-
02
8.567
e
-
02
ANFIS-QPSO-ADCEC
1.155
e
-
02
5.820
e
-
02Figure 6
The actual output and desired output by using the ANFIS-QPSO-ADCEC for Example5.


## 7. Conclusion

In this paper, we proposed a novel variant of the QPSO algorithm for training the parameters of ANFIS. This improved QPSO, called QPSO-ADCEC, employs an adaptive dynamical varying CE coefficient which makes the QPSO have better global search ability. The proposed QPSO-ADCEC is applied to estimation of the premise parameters of the ANFIS model, and the LSE approach is used for optimizing the consequent parameters of the fuzzy system.The effectiveness of the proposed ANFIS-QPSO-ADCEC method was verified by applying it to identification of nonlinear systems and predication of a chaotic system. The simulation results show that the proposed ANFIS-QPSO-ADCEC method had better performance than the ANFIS-PSO, ANFIS-QPSO, and the ANFIS trained with the gradient decent method due to the stronger global search ability of the QPSO-ADCEC algorithm.



---

*Source: 595639-2013-06-18.xml*