Identification of fractional-order systems via a switching differential evolution subject to noise perturbations

Physics Letters A 376 (2012) 3113–3120

Contents lists available at SciVerse ScienceDirect

Physics Letters A www.elsevier.com/locate/pla

Identification of fractional-order systems via a switching differential evolution subject to noise perturbations Wu Zhu a,1 , Jian-an Fang a,1 , Yang Tang b,c,d,∗ , Wenbing Zhang e , Yulong Xu a a

College of Information Science and Technology, Donghua University, Shanghai 201620, China Institute of Physics, Humboldt University, Berlin 12489, Germany c Potsdam Institute for Climate Impact Research, Potsdam 14415, Germany d Research Institute for Intelligent Control and System, Harbin Institute of Technology, Harbin 150006, China e Institute of Textiles and Clothing, The Hong Kong Polytechnic University, Hong Kong, China b

a r t i c l e

i n f o

Article history: Received 14 May 2012 Received in revised form 24 August 2012 Accepted 12 September 2012 Available online 3 October 2012 Communicated by A.R. Bishop Keywords: Fractional-order chaotic systems Differential evolution Population dynamics System identification

a b s t r a c t In this Letter, a differential evolution variant, called switching DE (SDE), has been employed to estimate the orders and parameters in incommensurate fractional-order chaotic systems. The proposed algorithm includes a switching population utilization strategy, where the population size is adjusted dynamically based on the solution-searching status. Thus, this adaptive control method realizes the identification of fractional-order Lorenz, Lü and Chen systems in both deterministic and stochastic environments, respectively. Numerical simulations are provided, where comparisons are made with five other Stateof-the-Art evolutionary algorithms (EAs) to verify the effectiveness of the proposed method. © 2012 Elsevier B.V. All rights reserved.

1. Introduction The theory of fractional calculus is a 300-year-old topic which can be traced back to Leibniz, Riemann, Liouville, Grünwald, and Letnikov [1]. For a long time, it did not attract much attention due to its inherent complexity and the fact that it does not have an acceptable geometrical or physical interpretation. However in recent years, more and more applications are found in different scientific fields, covering informatics, materials, physics and engineering [2–4]. It has been revealed that in interdisciplinary fields, various systems exhibit fractional dynamics [5,6]. For example, viscoelasticity, dielectric polarization, quantum evolution of complex system, fractional kinetics, and anomalous attenuation can be described by fractional differential equations [7,8]. On the other hand, the development of models based on chaotic systems has recently gained popularity in the investigation of dynamical systems because of its useful applications in many fields such as secure communication, data encryption, flow dynamics and biomedical engineering [9,10]. In the past years, chaos in integer-

*

Corresponding author at: Institute of Physics, Humboldt University, Berlin 12489, Germany. Tel.: +4933128820768; fax: +4930209399188. E-mail addresses: [email protected] (W. Zhu), [email protected] (Y. Tang). 1 These authors contributed equally to this work. 0375-9601/$ – see front matter © 2012 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.physleta.2012.09.042

order systems has been extensively studied. More recently, many researchers show growing interest in chaotic behavior of fractionalorder dynamical systems [11]. In Ref. [12], the authors found that chaotic behaviors exist in the fractional-order Rössler equation with orders less than 3, and hyperchaos exists in the fractionalorder Rössler hyperchaotic equation with order less than 4. Chaotic behavior of fractional-order Lorenz system was studied in Ref. [13]. In Ref. [14], it is demonstrated that the fractional-order Chua’s circuit of order as low as 2.7 can exhibit a chaotic attractor. In Ref. [10], chaotic behavior of fractional-order cellular neural networks was investigated, in which chaotic attractor is observed with the system orders as low as 2. As is well known, the complex dynamics of fractional-order chaotic systems has attracted increasing attention over recent years [3,8,11,15]. Among these literature, mainly two methods are used. One is optimization method [16,17] and the other is synchronization method [18–20]. Many synchronization methods are valid for fractional-order chaotic systems with known parameters. However, it is difficult to identify a given system from data when the physical systems are characterized by fractional-order differential equations with unknown orders or parameters. To the best of our knowledge, little results have been developed in the literature on parameter estimation of fractional-order chaotic systems based on chaotic synchronization. The difficulty lies in that the design of controller and the updating law of parameter identification

3114

W. Zhu et al. / Physics Letters A 376 (2012) 3113–3120

Fig. 1. Block diagram of the system identification process for a fractional-order chaotic system.

is a task with technique and sensitively depends on the considered systems. In search of an algorithm for parameters estimation of fractional-order chaotic systems, differential evolution (DE) algorithm [16,21] appears to be a competent candidate. Since its original definition by Storn and Price [22], the DE algorithm is perceived as a reliable and versatile population-based heuristic optimization technique, which exhibits remarkable performance in a wide variety of problems. In this Letter, we aim at adopting an improved differential evolution to identify the fractional-order Lorenz, Lü and Chen system in both deterministic and stochastic environments, respectively. The proposed algorithm is a kind of DE variants with a switching population utilization strategy (SDE). It considers the convergence speed and the computational cost simultaneously by nonperiodic partial increasing or declining individuals according to fitness diversities. Simulations have been done with the SDE-based evolutionary method to estimate the orders and parameters of incommensurate fractional-order chaotic systems. The performance is compared to that of three versions of DE, evolution strategy with covariance matrix adaptation (CMA-ES) and global and local real-coded genetic algorithm (GL-25). At the same time, we also analyze the basic properties of the switching approach. The rest of the Letter is organized as follows. In Section 2, objective function and the improved differential evolution variant are constructed. Numerical simulations are used to identify parameters and order of fractional-order chaotic systems in Section 3. In Section 4, a conclusion is drawn. 2. Parameters estimation of fractional-order chaotic systems using SDE 2.1. Problem formulation With the development of the regular calculus, fractional calculus is being known gradually. Fractional calculus is a generalization of the regular integer calculus. A positive order generalizes the common concepts of derivative and a negative order generalizes the common concepts of integral [1]. Unlike the numerical algorithm for solving an ordinary differential equation, the numerical simulation of a fractional differential equation is not so easy. In this Letter, we use the Caputo version and employ a predictor–corrector algorithm for fractional-order differential equations [23,24], which is the generalization of Adams–Bashforth– Moulton one. Consider the following n-dimensional fractionalorder chaotic systems:



D q Y (t ) = F Y (t ), t , θ



(1)

where Y (t ) = [ y 1 (t ), y 2 (t ), . . . , yn (t )] ∈ R indicates the n-dimensional state vector. θ = [θ1 , θ2 , . . . , θm ] T is a set of original parameters and q = [q1 , q2 , . . . , qn ] T ∈ R n denotes the fractional derivative orders. When q1 = q2 = · · · = qn , system (1) is called a commensurate ordered system, otherwise system (1) indicates an T

n

incommensurate ordered system. In this Letter, the fractionalorder q are also treated as additional parameters to be identified. Then, the generalized parameter vector of system (1) is Θ = [θ1 , θ2 , . . . , θm ; q1 , q2 , . . . , qn ] T . Suppose the structure of system (1) is known, thus the estimated system can be written as



D qˆ Yˆ (t ) = F Yˆ (t ), t , θˆ



(2)

where Yˆ (t ) = [Yˆ 1 (t ), Yˆ n (t ), . . . , Yˆ n (t )] T ∈ R n is the state vector of the estimated system. Let Θˆ = [θˆ1 , θˆ2 , . . . , θˆm ; qˆ 1 , qˆ 2 , . . . , qˆ n ] T be the estimated generalized parameter vector. To identify parameters, the following objective function is defined:

J (Θ) =

N 1 

N

e 2 (k) =

k =1

N 1 

N

2

Yˆ (k) − Y (k) ,

(3)

k =1

where Yˆ (k) and Y (k) are the actual and desired states of the system at time k, respectively; and e (k) = Yˆ (k) − Y (k) is the system’s error output; N denotes the length of data used for parameter estimation. The aim is to minimize the objective function J (Θ) by adjusting Θ , i.e.,

Θ ∗ = arg min J , Θ∈Ω

(4)

where Ω is searching space admitted for parameters and fractional-order. Clearly, the problem (4) is generally nonlinear and multi-modal, where the decision vector is Θ . The block diagram of the system identification process for a fractional-order chaotic system has been elaborated in Fig. 1. However, it is difficult to estimate the original parameters of the fractional-order chaotic systems using traditional optimization methods. Therefore, we employ an improved differential evolution (SDE) to find the true global minimum. 2.2. A switching differential evolution The stochastic system is iterated forward in time using a synchronous DE updating scheme. In our work, a mode-dependent population updating equation with Markovian switching parameters is introduced with the hope to keep track of the progress of individuals and further improve the search abilities. A detailed algorithm design of SDE can be found in Ref. [25]. If few trial vectors can outperform the corresponding parent in selection operation, particles may be clustered together and trapped into the local basin. In such a case, a sub-optimizer, called population-incremental strategy (PIS), is proposed to place several new individuals in appropriate areas to discover new possible solutions. However, in case most of individuals can spawn new promising offspings in the evolutionary process, it then signifies that redundant intermediate particles have been existed, and we introduce a sub-optimizer, namely population-cut strategy (PCS), to remove several poor particles with its entropy and ranking

W. Zhu et al. / Physics Letters A 376 (2012) 3113–3120

methods to reserve a place for a better reproduction. In SDE, the probability of selecting different sub-optimizers to improve the online solution-searching status is completely according to a nonhomogeneous Markov chain. For choosing required sub-optimizers adaptively, consider the following probability transition matrix in Eq. (5):



ϕa

ϕb ϕc

(5)

Here, M(1) = a, M(2) = b and M(3) = c stand for the population maintaining state, population increasing state and population decreasing state, respectively. The expectations of Markov chain Π are automatically updated by the search environment. It is worthwhile to mention that π22 and π33 are set to be 0. In this case, the incremental and cut operators will not be performed in consecutive generations. At each generation G, the conditional probability ϕb and ϕc is independently generated according to a Cauchy distribution with location parameters bootb and bootc , as well as standard deviation 0.1.

ϕb = Cauchy(bootb , 0.1),

(6)

ϕc = Cauchy(bootc , 0.1),

(7)

ϕa = 1 − ϕb − ϕc ,

(8)

and then truncated to be 1 if ϕb + ϕc > 1 or regenerated if ϕ  0. Denote boot as the location parameter of the Cauchy distribution that is initialized to be 0.1 and then updated at the end of each generation according to an exponential transformation. The flowchart of the DE algorithm (DE/rand/1/bin) with our proposed switching approach is elaborated in Fig. 3. In order to enhance the possibility of perturbation beyond the neighborhood to avoid local stagnation or premature convergence, the PIS is proposed to perturb the population and generate the necessary “fine” individuals.

 δ2 = a ·

(1 − λ(G − λ1 ))2 (1 − λ(G − λ1 ))2 + (1 + λ(G − λ1 ))2

 +b ,

(9)

where δ2 is the number of dimensions for perturbation, which is monotonically decreasing through the evolutionary search. Parameters a and b are the magnification coefficients. Parameters λ and G are here considered as the generation variables. During early stages of the optimization process, much more reproductions will be generated to spread out its particles within the decision space. Nevertheless, solutions of the population tend to concentrate in specific parts of the decision space during the later period of optimization process. On the other hand, the objective of the PCS is to prevent the excessive particles in current population and to reserve a space to accommodate a new potential particle. In SDE, the performance can easily be viewed by observing the fitness function value. In this case, two metrics will be adopted in objective space to keep track of exploration process. Then, the redundant (poor performing) particles from population should be removed stochastically.



τ = 1−

1 ranki + 1



Table 1 Problem illustration.

  × 1 − Hi (D) ,

(10)

where τ denotes an overall deletion indicator. The variable ranki and H i ( D ) indicate the rank metric and the entropy metric for individual i, respectively. It can be observed that from Eq. (10), for the individuals that have high rank values (i.e., away from the global best solution) or low entropy values (i.e., located in the crowded regions); these particles will have a higher probability of elimination.

ω(k)

Search space

N

0 and N (0, 1)

[0.5, 1.2] [30, 40] [−12, −2] [23, 33] [−2, 8]

100

q1 = 0.985, q2 = q3 = 0.99 σ = 30 ρ =0 α = 22.2

0 and N (0, 1)

[0.5, 1.2] [25, 35] [−5, 5] [17, 27] [−2, 8]

100

q1 = 0.985, q2 = q3 = 0.99 σ = 10 ρ = 28 α = −1

0 and N (0, 1)

[0.5, 1.2] [5, 15] [23, 33] [−6, 4] [−2, 8]

100

Instance

Θ

Chen system a=1

q1 = 0.985, q2 = q3 = 0.99 σ = 35 ρ = −7 α = 28 β =3



Π = 1 − ϕc 0 ϕc . 1 − ϕb ϕb 0

3115

Lü system a = 0.8

β=

Lorenz system a=0

β=

8.8 3

8 3

3. Simulation results 3.1. Experiments setup In 2002, Lü et al. [26] introduced a new chaotic system which can unify chaotic systems, such as, the Lorenz, Chen, Lü systems. The fractional-order differential equation of the unified system can be given by

⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩

d q1 x dt q1 d q2 x dt q2 d q3 x dt q3

= (25a + 10)( y − x), = (28 − 35a)x − xz + (29a − 1) y , = xy −

(11)

a +8 z, 3

where 0 < q i  1 (i = 1, 2, 3) is the order, a ∈ [0, 1]. Obviously, when a = 1, a = 0.8, a = 0, the unified chaotic system is the Chen, Lü and Lorenz system, respectively. The fractional-order unified chaotic system display chaotic attractors when (q1 , q2 , q3 ) = (0.985, 0.99, 0.99) [27]. Here, let us consider a more flexibility and universality replacement in the form of

⎧ σ = 25a + 10, ⎪ ⎪ ⎨ ρ = 28 − 35a,

α = 29a − 1, ⎪ ⎪ ⎩ β = a+3 8 ,

(12)

where σ , ρ , α and β are unknown parameters to be identified instead of a. With the addition of q i , the decision vector Θ has been generalized to be [q1 , q2 , q3 , σ , ρ , α , β] T . In this way, three typical system identification problems make up the test suite used for this comparative study, which are listed in Table 1. Θ specifies the generalized parameter vector; ω(k) is the system noise, which is independent of the input Y 0 ; N (0, 1) presents the white Gaussian noise (WGN) in zero-mean normal distribution with variance 1; Search space is the predefined boundary constraints; N denotes the data length used in calculating the mean square error (MSE). For each algorithm and each test function, 30 independent runs are conducted with D ∗ 10 000 FES as the termination criterion, while D = 7 defines the number of decision space. Traditionally, “generation” is a natural form of computational cost for statistical comparison. However, the population may not be the same in different algorithms. The algorithm with a larger population may obtain a better performance together with much more function evaluations in every generation. Thus, in this Letter, the function evaluations (FES) are conducted to represent its computational cost for algorithm comparison.

3116

W. Zhu et al. / Physics Letters A 376 (2012) 3113–3120

Table 2 EA algorithms for comparison. Algorithm

Parameters

Reference

GL-25 CMA-ES jDE SaDE CoDE SDE

P G = 25%, N FL = 5 σ = 0.25 τ 1 = τ 2 = 0.1 F : rand N (0.5, 0.3), L P = 50 automatically chosen p 1 = p 2 = 1, ζ = 4

[35] [36] [37] [28] [29] this Letter

In all simulations, the population size (SP) of the most EAs is 100 with the exception of SaDE and CoDE. As suggested in Ref. [28], the population size of SaDE is chosen to be 50. And SP of CoDE is set as 30 in Ref. [29]. Five existing EA algorithms are shown in Table 2 in detail. GL-25 is a hybrid real-coded genetic algorithm which combines the global and local search. CMA-ES represents the State-of-the-Art of evolution strategies and it is a referent in the continuous optimization field. jDE is a standard DE with adapted parameter setting. SaDE delivers a mutation strategy pool where strategy is self-adapted based on their previous performance. CoDE is an adaptive DE, where three mutation strategies are coupled with three parameter settings randomly to generate offsprings. The parameters for these EAs are provided in Table 2. 3.2. Comparisons on the solution accuracy In this section, an overall comparison of the performance is provided between the SDE variant and other five State-of-the-Art EAs (i.e., GL-25, CMA-ES, jDE, SaDE and CoDE). We evaluate the performance of the six heuristic algorithms over three typical nonlinear uncertain discrete-time systems. Fig. 4 illustrates the objective function value versus number of evaluations averaged over 30 random runs for the six algorithms. Table 3 reports the experimental results of three systems in deterministic environment, averaged over 30 independent runs with D ∗ 10 000 FES. The best result

among those EAs is indicated by Boldface in the table. Fig. 2 also shows instance of evolution of the parameters of three systems for SDE. From Table 3 and Fig. 4, the SDE provides the best performance on a = 0.8 and a = 0, then ranks the second on a = 1. jDE gives the best performance on a = 1. The results show that CMA-ES and jDE have good ability of convergence speed. To be specific, on the Chen system (a = 1), jDE delivers good accuracy and the highest convergence rate, while SDE outperforms other four methods. On the Lü system (a = 0.8), SDE performs significantly better than five others. The outstanding performance of SDE is due to its dynamic PCS, which leads to very fast convergence. On the Lorenz system (a = 0), SDE and CMA-ES achieve similar performance and outperform other four methods. For a thorough comparison, the t-test has also been carried out in this Letter. Table 3 presents the score on every function of this two-tailed test with a significance level of 0.05 between the SDE variant and other heuristic algorithms. Marks “+ (Better)”, “≈ (Same)” and “− (Worse)” denote that the performance of SDE is significantly better than, almost the same as, and significantly worse than the compared algorithm on fitness values in 30 runs, respectively. As confirmed in t-test, the SDE in general offers more improved performance than the other five State-of-the-Art EAs. As mentioned above, the SDE has shown a very competitive performance in the three fractional-order chaotic systems. In practical engineering, noise exists universally in nature [30,31]. Therefore, in the past few decades researchers have witnessed significant progress on filtering and control for linear/nonlinear systems with various types of noises, among which Gaussian noise is one of the most general signals that has been widely studied [32,33]. Here, we further evaluate the proposed framework on the three expanded stochastic systems, where a zero-mean Gaussian whitenoise is added. The maximum number of FES is set to be D ∗ 10 000 in all runs. Table 4 summarizes the experimental results. Fig. 4 illustrates the objective function value versus number of evaluations averaged over 30 random runs for the six algorithms.

Table 3 Experimental results in deterministic environment, averaged over 30 independent runs with D ∗ 10 000 FES. Instance

a=1 Mean error ± std. dev.

a = 0.8 Mean error ± std. dev.

a=0 Mean error ± std. dev.

GL-25 CMA-ES jDE SaDE CoDE SDE

3.71E−03±2.32E−03+ 1.51E−02±4.73E−02+ 1.14E−05±2.14E−05− 5.79E−04±7.77E−04+ 1.56E−04±2.85E−04+ 4.72E−05±1.75E−04

5.43E−04±9.13E−04+ 2.12E−04±8.03E−04+ 1.30E−06±1.02E−06+ 5.26E−04±1.01E−03+ 4.32E−04±8.83E−04+ 6.13E−07±1.09E−13

2.86E−05±3.99E−05+ 2.58E−26±2.44E−26+ 2.38E−07±4.89E−07+ 3.10E−07±7.64E−07+ 3.07E−12±1.68E−11+ 1.11E−27±2.89E−28

Fig. 2. Instance of evolving process of parameters identification of fractional-order Lorenz system in deterministic environment.

W. Zhu et al. / Physics Letters A 376 (2012) 3113–3120

3117

Fig. 3. Flowchart of the DE algorithm with the switching population utilization strategy.

From Table 4 and Fig. 4, the SDE provides the best performance on a = 1, a = 0.8 and a = 0. jDE offers the same performance on a = 0.8. The results show that SaDE and jDE have a good ability of convergence speed. Hence, SDE exhibits the highest performance in noise-expanded identification problems, which can efficiently adjust the population structure and guide the evolution process toward more promising solutions. The t-test is also illustrated in Table 4. In fact, SDE performs better than GL-25, CMA-ES, jDE, SaDE and CoDE on 3, 3, 2, 3 and 3 out of 3 test functions, respectively. Thus, SDE is better than other five competitors in fractional-order chaotic system identification problems. Comparing the results and the cost function graphs, among these six EAs, the CMA-ES can converge to the best solution found so far very quickly though it is easy to stuck in the local optima. The SaDE has a good global search ability and slow convergence speed. The jDE has a good search capability on noise-expanded filtering problems. The SDE has a

good local search ability and global search ability at the same time. 3.3. Comparison on convergence rate and successful percentage The convergence rate for achieving the global optimum is another key point for testing the algorithm performance. Note that in solving real-world optimization problems, the “function evaluation” overwhelms the algorithm overhead. Hence, the computation times of these algorithms are not provided here. Table 5 shows that SDE needs the least FES to achieve the acceptable solution on a = 0.8 in deterministic environment and a = 1 in stochastic environment, which reveals that SDE has a higher convergence rate than other algorithms. Though CMA-ES or GL-25 might outperform SDE on the other functions, CMA-ES and GL-25 have much worse successful ratio and accuracy than SDE on the problems for comparison. In addition, SDE can achieve accepted value with a

3118

W. Zhu et al. / Physics Letters A 376 (2012) 3113–3120

Fig. 4. Cost function value versus number of evaluations averaged over 30 random runs for the six algorithms. In deterministic environment: (a) a = 1, (b) a = 0.8, (c) a = 0. In stochastic environment: (d) a = 1, (e) a = 0.8, (f) a = 0.

good convergence speed and accuracy on most of the problems, as shown in Tables 3–4 and Fig. 4. Tables 3–4 also show that SDE yields a highest percentage for achieving acceptable solutions in 30 runs. According to the no free lunch theorem [34], any elevated performance over one class of problems is offset by performance over another class. Hence, one algorithm cannot perform better on convergence speed and accuracy than the others on every optimization problem. In summary, the SDE performs best on parameter identification problems with or without noise. Owing to the switching technique,

the SDE processes capabilities of fast convergence speed, highest successful ratio and the best search accuracy among these EAs. 3.4. Performance of the switching approach In this section, the switching approach is used to test the search performance of SDE. In all the experiments, threshold ζ is adjusted in the following. Moreover, the parameter p 2 is also considered, which denotes the number of potential candidates for perturbation.

W. Zhu et al. / Physics Letters A 376 (2012) 3113–3120

3119

Table 4 Experimental results in stochastic environment, averaged over 30 independent runs with D ∗ 10 000 FES. Instance

a=1 Mean error ± std. dev.

a = 0.8 Mean error ± std. dev.

a=0 Mean error ± std. dev.

GL-25 CMA-ES jDE SaDE CoDE SDE

4.69E−01±4.17E−02+ 9.91E−01±6.84E−01+ 4.35E−01±3.18E−02+ 4.86E−01±2.79E−02+ 5.08E−01±4.13E−02+ 4.19E−01±3.40E−02

4.47E−01±3.58E−02+ 1.16E+00±1.79E+00+ 4.30E−01±4.57E−02≈ 4.93E−01±4.94E−02+ 5.21E−01±4.69E−02+ 4.30E−01±3.41E−02

4.24E−01±3.71E−02+ 1.05E+00±9.88E−01+ 4.43E−01±2.89E−02+ 4.14E−01±2.65E−02+ 4.51E−01±4.83E−02+ 4.06E−01±3.15E−02

Table 5 Convergence speed and algorithm reliability comparisons, averaged over 30 independent runs with D ∗ 10 000 FES. Instance

Deterministic environment

Stochastic environment

a=1

a = 0.8

a=0

a=1

a = 0.8

a=0

–

–

24 726 86.7

22 470 70

26 571 90

–

–

GL-25

Mean FEs Right percentage (%)

38 407 6 .7

CMA-ES

Mean FEs Right percentage (%)

5337 80

48 420 96.7

14 763 100

–

Mean FEs Right percentage (%)

28 190 100

61 414 23.3

–

27 242 86.7

31 355 96.7

26 372 73.3

Mean FEs Right percentage (%)

37 283 70

– 0

51 350 3.3

35 157 70

36 779 40

21 994 96.7

Mean FEs Right percentage (%)

33 283 96.7

– 0

69 060 3 .3

40 995 40

35 610 20

29 855 60

Mean FEs Right percentage (%)

24 103 96.7

23 366 100

60 241 100

23 191 100

32 437 100

31 562 100

jDE

SaDE

CoDE

SDE

0

0

0

0

0

0

Table 6 Effects of ζ and p 2 on search accuracy of SDE. (The best value is in bold font.) Instance Theoretical condition a=1 a = 0.8 a=0 Stochastic condition a=1 a = 0.8 a=0

ζ = 4, p 2 = 5 Mean error ± std. dev.

ζ = 2, p 2 = 1 Mean error ± std. dev.

ζ = 4, p 2 = 1 Mean error ± std. dev.

2.71E−05±1.41E−04≈ 1.06E−04±5.78E−04+ 1.62E−16±5.41E−16+

3.95E−04±1.25E−03+ 6.14E−07±2.56E−10≈ 4.23E−25±2.31E−24+

4.72E−05±1.75E−04 6.13E−07±1.09E−13 1.11E−27±2.89E−28

4.39E−01±3.09E−02+ 4.26E−01±3.04E−02≈ 4.31E−01±2.96E−02+

4.26E−01±3.78E−02≈ 4.30E−01±3.17E−02≈ 4.13E−01±3.62E−02+

4.19E−01±3.40E−02 4.30E−01±3.41E−02 4.06E−01±3.15E−02

In this Letter, ζ indicates the trigger thresholds, which is used to control the sensitivity of the dynamic SDE. While ζ is set as one, the population size will be adjusted in each generation. Setting a higher ζ value will result in a lower sensitivity of the SDE, while a lower ζ value will issue in a higher efficiency of the population adjustment. Notice that the parameter ζ should set to be larger than one. Failure to do this will result in an instant elimination of a newborn individual with poor performance, which may provide some degree of diversity preservation. On the other hand, coefficient p 2 also influences the perturbation process substantially. Table 6 shows the comparisons between SDE with other three parameter settings of SDE over three parameter identification problems with or without noise perturbation. It indicates that SDE is not sensitive to the adjustment of parameters. In order to make a balance of the search accuracy and robustness, ζ = 4 and p 2 = 1 are used as a representative parameter setting in our Letter. This setting will prevent the instant elimination of a newborn individual and keep the dynamic SDE high sensitivity.

ture of SDE is the switching population utilization strategy which improves the quality of the population and decreases the number of calculations by nonperiodic partial increasing or declining individuals. The results are shown in comparison with five other existing methods. The results obtained by our approach are better than other EAs especially in stochastic environment. We believe that such a simple adaptive scheme will be much beneficial for the applications of fractional order systems identification. Acknowledgements This work was supported by the Fundamental Research Funds for the Central Universities, the Alexander von Humboldt Foundation of Germany, the National Natural Science Foundation of PR China (61203235), the China Postdoctoral Science Foundation (2011M501040) and the Key Creative Project of Shanghai Education Community (13ZZ050). References

4. Conclusion This Letter presents a switching differential evolution algorithm for identifying the orders and parameters of incommensurate fractional-order Lorenz, Lü and Chen systems. The main fea-

[1] [2] [3] [4]

I. Podlubny, Fractional Differential Equations, Academic Press, New York, 1999. G. Wu, S. Zhang, Phys. Lett. A 375 (2011) 3659. Y. Tang, Z. Wang, J. Fang, Chaos 19 (2009) 013112. R. Hilfer, Applications of Fractional Calculus in Physics, World Scientific, Singapore, 2000.

3120

W. Zhu et al. / Physics Letters A 376 (2012) 3113–3120

[5] G. Zaslavsky, Phys. Rep. 371 (2002) 461. [6] S. Hempel, A. Koseska, J. Kurths, Z. Nikoloski, Phys. Rev. Lett. 107 (2011) 054101. [7] J. He, S. Elagan, Z. Li, Phys. Lett. A 376 (2012) 257. [8] L. Yuan, Q. Yang, Commun. Nonlinear. Sci. Numer. Simulat. 17 (2012) 305. [9] H. Eleuch, A. Prasad, Phys. Lett. A 376 (2012) 1970. [10] I. Grigorenko, E. Grigorenko, Phys. Rev. Lett. 91 (2003) 034101. [11] Y. Tang, X. Zhang, C. Hua, L. Li, Y. Yang, Phys. Lett. A 376 (2012) 457. [12] C. Li, G. Chen, Physica A 341 (2004) 55. [13] P. Arena, L. Fortuna, D. Porto, Phys. Rev. E 61 (2000) 776. [14] T. Hartley, C. Lorenzo, H. Qammer, IEEE T. Circuits-I 42 (1995) 485. [15] Y. Zheng, Y. Nian, D. Wang, Phys. Lett. A 375 (2010) 125. [16] S. Das, P.N. Suganthan, IEEE T. Evolut. Comput. 15 (2011) 4. [17] A. Arenas, A. Díaz-Guilera, J. Kurths, Y. Moreno, C. Zhou, Phys. Rep. 469 (2008) 93. [18] Y. Tang, R. Qiu, J. Fang, Q. Miao, M. Xia, Phys. Lett. A 372 (2008) 4425. [19] Y. Tang, Z. Wang, W.K. Wong, J. Kurths, J. Fang, Chaos 21 (2011) 025114. [20] U. Vincent, R. Guo, Phys. Lett. A 375 (2011) 2322. [21] W. Zhu, J. Fang, Y. Tang, W. Zhang, W. Du, PLoS ONE 7 (2012) e40549. [22] R. Storn, K. Price, J. Global. Optim. 11 (1997) 341.

[23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37]

K. Diethelm, N. Ford, A. Freed, Nonlinear Dynam. 29 (2002) 3. D. Kunsezov, A. Bulagc, G. Dang, Phys. Rev. Lett. 82 (1999) 1136. W. Zhu, J. Fang, Y. Tang, W. Zhang, PLoS ONE, submitted for publication. J. Lu, G. Chen, D. Cheng, S. Celikovsky, Int. J. Bifurcat. Chaos 12 (2002) 2917. L. Sheu, H. Chen, J. Chen, L. Tam, W. Chen, K. Lin, Y. Kang, Chaos Soliton Fract. 36 (2008) 98. A.K. Qin, V.L. Huang, P.N. Suganthan, IEEE T. Evolut. Comput. 13 (2009) 398. Y. Wang, Z. Cai, Q. Zhang, IEEE T. Evolut. Comput. 15 (2011) 55. Y. Tang, W. Zou, J. Lu, J. Kurths, Phys. Rev. E 85 (2012) 046207. K. Somiya, A.G. Gurkovsky, D. Heinert, S. Hild, R. Nawrodt, S. Vyatchanin, Phys. Lett. A 375 (2011) 1363. H. Peng, L. Li, Y. Yang, F. Sun, Phys. Rev. E 83 (2011) 036202. L. Lafuerza, P. Colet, R. Toral, Phys. Rev. Lett. 105 (2010) 084101. D.H. Wolpert, W.G. Macready, IEEE T. Evolut. Comput. 1 (1997) 67. C. Garcia-Martinez, M. Lozano, F. Herrera, D. Molina, A.M. Sanchez, Eur. J. Oper. Res. 185 (2008) 1088. N. Hansen, A. Ostermeier, Evol. Comput. 9 (2001) 159. J. Brest, S. Greiner, B. Boscovic, M. Mernik, V. Zumer, IEEE T. Evolut. Comput. 10 (2006) 646.