Completely Derandomized Self-Adaptation in Evolution
Strategies
Nikolaus Hansen and Andreas Ostermeier
In Evolutionary Computation, 9(2), pp. 159-195 (2001)
Errata• Section 3, footnote 9: ”We use the expectation of nk ln(...)” must be ”We use the expectationof 21 nk ln(...)”.
• (same error) Section 3, same page, description Performance: ”...where ln(...)...” must be”...where 21 ln(...)...”.
Additional Comments• Section 5.1, Table 1: The default parameter setting for wi makes sense only, if µ ≤ λ/2.
Otherwise weights become negative which is unjustified in our context and was excluded bydefinition. You may choose wi = ln(max(λ/2, µ) + 1/2)) − ln(i) instead. This setting isapplicable to any µ ≤ λ.
• Section 5.1, discussion of parameter cc : ”..., we suspect cσ ≤ cc ≤ 1 to be a sensible choicefor cc .” In fact, up to now, we found no particular evidence against choosing cc smaller thancσ . In any case it is necessary to choose cc > ccov .
• Section 6, Table 3: the input vector argument to functions 2-12 on the left side of the equations
should rather be x instead of y.

Completely Derandomized Self-Adaptation in
Evolution Strategies
Nikolaus Hansen
hansen@bionik.tu-berlin.deTechnische Universität Berlin, Fachgebiet für Bionik, Sekr. ACK 1, Ackerstr. 71–76,
13355 Berlin, Germany
Andreas Ostermeier

ostermeier@bionik.tu-berlin.deTechnische Universität Berlin, Fachgebiet für Bionik, Sekr. ACK 1, Ackerstr. 71–76,
13355 Berlin, Germany

Abstract
This paper puts forward two useful methods for self-adaptation of the mutation distribution – the concepts of derandomization and cumulation. Principle shortcomings ofthe concept of mutative strategy parameter control and two levels of derandomization
are reviewed. Basic demands on the self-adaptation of arbitrary (normal) mutation distributions are developed. Applying arbitrary, normal mutation distributions is equivalent to applying a general, linear problem encoding.
The underlying objective of mutative strategy parameter control is roughly to favorpreviously selected mutation steps in the future. If this objective is pursued rigorously, a completely derandomized self-adaptation scheme results, which adapts arbitrary normal mutation distributions. This scheme, called covariance matrix adaptation(CMA), meets the previously stated demands. It can still be considerably improved bycumulation – utilizing an evolution path rather than single search steps.
Simulations on various test functions reveal local and global search properties ofthe evolution strategy with and without covariance matrix adaptation. Their performances are comparable only on perfectly scaled functions. On badly scaled, nonseparable functions usually a speed up factor of several orders of magnitude is observed. On moderately mis-scaled functions a speed up factor of three to ten can be
expected.
Keywords
Evolution strategy, self-adaptation, strategy parameter control, step size control, derandomization, derandomized self-adaptation, covariance matrix adaptation, evolution path, cumulation, cumulative path length control.

1Introduction

The evolution strategy (ES) is a stochastic search algorithm that addresses the following
search problem: Minimize a non-linear objective function that is a mapping from search
spaceto . Search steps are taken by stochastic variation, so-called mutation, of(recombinations of) points found so far. The best out of a number of new search points
are selected to continue. The mutation is usually carried out by adding a realizationof a normally distributed random vector. It is easy to imagine that the parametersof the normal distribution play an essential role for the performance1 of the search 

1 Performance, as used in this paper, always refers to the (expected) number of required objective function
evaluations to reach a certain function value.

c 2001 by the Massachusetts Institute of Technology

Evolutionary Computation 9(2): 159-195

N. Hansen and A. Ostermeier

Figure 1: One- lines of equal probability density of two normal distributions respectively. Left: one free parameter (circles). Middle:  free parameters (axis-parallel ellipsoids). Right: 
free parameters (arbitrarily oriented ellipsoids).
algorithm. This paper is specifically concerned with the adjustment of the parametersof the normal mutation distribution.
Among others, the parameters that parameterize the mutation distribution are
called strategy parameters, in contrast to object parameters that define points in search
space. Usually, no particularly detailed knowledge about the suitable choice of strategy parameters is available. With respect to the mutation distribution, there is typically
only a small width of strategy parameter settings where substantial search progress
can be observed (Rechenberg, 1973). Good parameter settings differ remarkably fromproblem to problem. Even worse, they usually change during the search process (possibly by several orders of magnitude). For this reason, self-adaptation of the mutation
distribution that dynamically adapts strategy parameters during the search process isan essential feature in ESs.
We briefly review three consecutive steps of adapting normal mutation distributions in ESs.
1. The normal distribution is chosen to be isotropic. Surfaces of equal probability
density are circles (Figure 1, left), or (hyper-)spheres if  . Overall variance ofthe distribution –– in other words the (global) step size or the expected step length ––is the only free strategy parameter.
2. The concept of global step size can be generalized.2 Each coordinate axis is assigned a different variance (Figure 1, middle) –– often referred to as individual step
sizes. There are  free strategy parameters. The disadvantage of this concept is thedependency on the coordinate system. Invariance against rotation of the search
space is lost. Why invariance is an important feature of an ES is discussed in Section 6.
3. A further generalization dynamically adapts the orthogonal coordinate system,
where each coordinate axis is assigned a different variance (Figure 1, right). Any
normal distribution (with zero mean) can be produced. This concept results in
2 There is more than one sensible generalization. For example, it is possible to provide one arbitrarilyoriented axis with a different variance. Then  parameters have to be adapted. Such an adaptation canbe formulated independent of the coordinate system.

2

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation    free strategy parameters. If an adaptation mechanism is suitably formulated, the invariance against rotations of search space is restored.
The adaptation of strategy parameters in ESs typically takes place in the concept ofmutative strategy parameter control (MSC). Strategy parameters are mutated, and a new
search point is generated by means of this mutated  strategy parameter setting.
We exemplify the concept formulating an   -ES with (purely) mutative control
of one global step size.  and  
are the object parameter vector and
step size3 of the parent at generation  . The mutation step from generation  to  
reads for each offspring 






 
where:# , for &






$#
 "!    %    

(1)
(2)

'(independent
realizations of a random number with#zero mean. Typically,  is normally distributed with standard
+#  deviation
and Schwefel, 1993). We usually prefer to choose *')$#  (Bäck
-,  .(Rechenberg, 1994).
*
0/1,  2

% 4365 8 7 : 9 ; , for < '( independent realizations of a 87 :9  -normally
9

distributed random vector, where  is the unity matrix. That is, components of%  are independent and identically ,   -normally distributed.
After mutation steps are carried out, the best offspring (with respect to the object parameter vector    ) is selected to start the next generation step. Equation (1)#facilitates the mutation on the strategy parameter level. The standard deviation ofrepresents the mutation strength on the strategy parameter level.
This adaptation concept was introduced by Rechenberg (1973) for global and individual step sizes. Schwefel (1981) expanded the mutative approach to the adaptationof arbitrary normal mutation distributions, which is discussed in more detail in Section
3.1. Ostermeier et al. (1994b) introduced a first level of derandomization into strategy
parameter control that facilitates an individual step size adaptation in constantly small
populations. (If the standard mutative strategy parameter control is used to adapt individual step sizes, based on our experience, the population size has to scale linearely
with the problem dimension  .)In this paper, a second level of derandomization is put forward: The original objective of mutative strategy parameter control –– that is to favor strategy parameter settings
that produce selected steps with high probability (again) –– is explicitely realized. Complete derandomization, applied to the adaptation of arbitrary normal mutation distributions, leads almost inevitably to the covariance matrix adaptation (CMA) described in
Section 3.2.
The paper is organized as follows. In Section 2, we discuss the basic shortcomingsof the concept of mutative strategy parameter control and review the derandomized
approach to strategy parameter control in detail. In Section 3, the different concepts areinvestigated with respect to the adaptation of arbitrary normal mutation distributions.
3 In this paper, step size always refers to = . Step size = is the component-wise standard deviationofthevector =>@?BA"C8DEGFIHKJ . For this vector, =L can be interpreted as overall variance andMNPO =>@random?BA"C8D'E ORQS =T  is the expected step length.

Evolutionary Computation Volume 9, Number 2

3

N. Hansen and A. Ostermeier

General demands on such an adaptation mechanism are developed, and the CMA-ES
is introduced. In Section 4, the objective of strategy parameter control is expanded
using search paths (evolution paths) rather than single search steps as adaptation criterion. This concept is implemented by means of so-called cumulation. In Section 5,
we formulate the CMA-ES algorithm, an evolution strategy that adapts arbitrary, normal mutation distributions within a completely derandomized adaptation scheme withcumulation. Section 6 discusses the test functions used in Section 7, where various simulation results are presented. In Section 8, a conclusion is given. Appendix A provides
a MATLAB implementation of the CMA-ES.

2

Derandomization of Mutative Strategy Parameter Control

In the concept of mutative strategy parameter control (MSC), the selection probability of a strategy parameter setting is the probability to produce (with these strategy
parameters) an object parameter setting that will be selected. Selection probability ofa strategy parameter setting can also be identified with its fitness. Consequently, we
assume the following aim behind the concept of MSC:4 Find the strategy parameterswith the highest selection probability or, in other words, raise the probability of mutation steps that produced selected individuals.5 This implicitly assumes that strategy
parameters that produced selected individuals before are suitable parameters in the(near) future. One idea of derandomization is to increase the probability of producing
previously selected mutation steps again in a more direct way than MSC.
Before the concept of derandomization is discussed in detail, some important
points concerning the concept of MSC are reviewed.
Selection of the strategy parameter setting is indirect. The selection process operates on the object parameter adjustment. Comparing two different strategy parameter settings, the better one has (only) a higher probability to be selected –– due
to object parameter realization. Differences between these selection probabilities
can be quite small, that is, the selection process on the strategy parameter level ishighly disturbed. One idea of derandomization is to reduce or even eliminate this
disturbance.
The mutation on the strategy parameter level, as with any mutation, produces different individuals that undergo selection. Mutation strength on the strategy parameter level must ensure a significant selection difference between the individuals.
In our view, this is the primary task of the mutation operator.
Mutation strength on the strategy parameter level is usually kept constantthroughout the search process. Therefore, the mutation operator (on the strategy
parameter level) must facilitate an effective mutation strength that is virtually independent of the actual position in strategy parameter space. This can be complicated, because the best distance measure may not be obvious, and the positionindependent formulation of a mutation operator can be difficult (see Section 3.1).
4 We assume there is one best strategy parameter setting at each time step. Alternatively, to apply different
parameter settings at the same time means to challenge the idea of a normal search distribution. We foundthis alternative to be disadvantageous.
5 Maximizing the selection probability is, in general, not identical with maximizing a progress rate. This is
one reason for the often observed phenomenon that the global step size is adapted to be too small by MSC.
This problem is addressed in the path length control (Equations (16) and (17)) by the so-called cumulation(Section 4).

4

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation

The possible (and the realized) change rate of the strategy parameters between
two generations is an important factor. It gives an upper limit for the adaptation
speed. When adapting  or even  strategy parameters, where  is the problem
dimension, adaptation speed becomes an important factor for the performance of
the ES. If only one global step size is adapted, the performance is limited by thepossible adaptation speed only on a linear objective function (see discussion ofparameter from (3) below).
On the one hand, it seems desirable to realize change rates that are as large as
possible –– achieving a fast change and consequently a short adaptation time. On
the other hand, there is an upper bound to the realized change rates due to thefinite amount of selection information. Greater changes cannot rely on valid information and lead to stochastic behavior. (This holds for any adaptation mechanism.) As a simple consequence, the change rate of a single strategy parametermust decrease with an increasing number of strategy parameters to be adapted,
assuming a certain constant selection scheme.
In MSC, the change rate of strategy parameters obviously depends (more or lessdirectly) on the mutation strength on the strategy parameter level. Based on thisobservation, considerable theoretical efforts were made to calculate the optimal
mutation strength for the global step size (Schwefel, 1995; Beyer, 1996b). But, ingeneral, the conflict between an optimal change rate versus a significant selection
difference (see above) cannot be resolved by choosing an ambiguous compromisefor the mutation strength on the strategy parameter level (Ostermeier et al., 1994a).
The mutation strength that achieves an optimal change rate is usually smaller thanthe mutation strength that achieves a suitable selection difference. The discrepancy increases with increasing problem dimension and with increasing number ofstrategy parameters to be adapted. One idea of derandomization is to explicitly
unlink the change rate from the mutation strength resolving this discrepancy.
Parent number and the recombination procedure have a great influence on thepossible change rate of the strategy parameters between (two) generations. Assuming a certain mutation strength on the strategy parameter level, the possible
change rate can be tuned downwards by increasing . This is most obvious forintermediate multi-recombination: The mean change of recombined individuals is approximately )times smaller than the mean change of a single individual. Within the concept of MSC, choosing and an appropriate recombination
mechanism is the only way to tune the change rate independently of the mutation
strength (downwards). Therefore, it is not a surprising observation that a successful strategy parameter adaptation in the concept of MSC strongly depends on asuitable choice of : In our experience has to scale linearly with the numberof strategy parameters to be adapted. One objective of derandomization is to facilitate a reliable and fast adaptation of strategy parameters independent of the
population size, even in small populations.














We start our discussion of derandomization from the ES with mutative strategy
parameter control formulated in Equations (1) and (2). For the sake of simplicity, we
still consider adaptation of one global step size only. (The technique of derandomization becomes especially important if a large number of strategy parameters has to beadapted.) The first level of derandomization (Ostermeier et al., 1994a) facilitates independent control of mutation strength and change rate of strategy parameters. This canEvolutionary Computation Volume 9, Number 2

5

N. Hansen and A. Ostermeier'(be achieved by slight changes in the formulation of the mutation step. For 


 









!



 



#



 "!

,(3)

+#  %  


(4)

where
is the damping parameter, and the symbols from Equations (1) and (2)in Section 1 are used. This can still be regarded as a mutative approach to strategy
parameter control: If 
, Equations (3) and (4) are identical to (1) and (2), because

8# 
 "!
 in Equation (4) equals  in Equation (2). Enlarging the damping para
meter reduces the change rate between  and  leaving the selection# relevant
mutation strength on the strategy parameter level (the standard deviation of  in Equation (4)) unchanged. Instead of choosing the standard deviation to be )  and ,
as in the purely mutative approach, it is more sensible to choose the standard deviationand , facilitating the same rate of change with a larger
 , assuming 
mutation strength. Large values for scale down fluctuations that occur due to stochas
tic effects, but decrease the possible change rate. For standard deviation
within theinterval , neither the fluctuations for small nor the limited change rate.
for large decisively worsen the performance on the sphere objective function
But if  or even  strategy parameters have to be adapted, the trade-off between large
fluctuations and a large change rate come the fore. Then the change rate must be tuned
carefully.
As Rechenberg (1994) and Beyer (1998) pointed out, this concept –– mutating big
and inheriting small –– can easily be applied to the object parameter mutation as well.
It will be useful in the case of distorted selection on the object parameter level and is
implicitly implemented by intermediate multi-recombination.
In general, the first level of derandomization yields the following effects:



 





  

 #Mutation strength (in our example, standard deviation of  ) can be chosen comparatively large in order to reduce the disturbance of strategy parameter selection.
The ratio between change rate and mutation strength can be tuned downwards (inour example, by means of ) according to the problem dimension and the numberof strategy parameters to be adapted. This scales down stochastic fluctuations of
the strategy parameters.

and can be chosen independent of the adaptation mechanism. They, in particular, become independent of the number of strategy parameters to be adapted.
The second level of derandomization completely removes the selection disturbance on the strategy parameter level. To achieve this, mutationon the strategy para#
meter level by means of sampling the random number  in Equation (3) is omitted.
Instead, the realized step length % is used to change . This leads to

 


 









"!

 



 %   / 5  7 (9  

(5)

 % 
#

(6)

 %   / 5 87 (9   in Equation (5) replaces in Equation (3). Note that  5 +7 :9  
is the expectation of  %   under random selection.6 Equation (5) is not intended to
6 Random selection occurs if the objective function returns random numbers independent of

6!.

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation 

be used in practice, because the selection relevance of %  vanishes for increasing  .7
One can interpret Equation (5) as a special mutation on the strategy parameter levelor as an adaptation procedure for (which must only be done for selected individuals). Now, all stochastic variation of the object parameters –– namely originated by the
random vector realization %  in Equation (6) –– is used for strategy parameter adaptation in Equation (5). In other words, the actually realized mutation step on the object
parameter level is used for strategy parameter adjustment. If %  is smaller than expected,  is decreased. If %  is larger than expected   is increased. Selecting asmall/large mutation step directly leads to reduction/enlargement of the step size.
The disadvantage of Equation (5) compared to Equation (1) is that it cannot be
expanded easily for the adaptation of other strategy parameters. The expansion to theadaptation of all distribution parameters is the subject of this paper.
In general, the completely derandomized adaptation follows the following principles:

 

 

1. The mutation distribution is changed, so that the probability to produce the selected mutation step (again) is increased.
2. There is an explicit control of the change rate of strategy parameters (in our example, by means of ).
3. Under random selection,
strategy parameters are stationary. In our example, exequals   .
pectation of  
In the next section, we discuss the adaptation of all variances and covariances of themutation distribution.

3

Adapting Arbitrary Normal Mutation Distributions

We give two motivations for the adaptation of arbitrary normal mutation distributions
with zero mean.
The suitable encoding of a given problem is crucial for an efficient search withan evolutionary algorithm.8 Desirable is an adaptive encoding mechanism. To
restrict such an adaptive encoding/decoding
to linear transformations seems reasonable: Even in the linear case,    parameters have to be adapted. It is easy
to show that in the case of additive mutation, linear transformation of the search
space (and search points accordingly) is equivalent to  linear
7 (9 transformation of themutation distribution. Linear transformation of the -normal mutation distribution always yields a normal distribution with zero mean, while any normal
distribution with zero mean can be produced by a suitable linear transformation.
This yields equivalence between an adaptive general linear encoding/decoding
and the adaptation of all distribution parameters in the covariance matrix.
We assume the objective function to be non-linear and significantly nonseparable –– otherwise the search problem is comparatively easy, because it canbe solved by solving  one-dimensional problems. Adaptation to non-separable
objective functions is not possible if only individual step sizes are adapted. This
7 For an efficient, completely derandomized global step size adaptation, the evolution path  ?  E     ? E  , where  ?  !E#"$%"&  ?  :E , must be used instead of   in Equation (5).

This method is referred to as cumulative path length control and implemented with Equations (16) and (17).
8 In biology, this is equivalent to a suitable genotype-phenotype transformation.
Evolutionary Computation Volume 9, Number 2

7

N. Hansen and A. Ostermeierdemands a more general mutation distribution and suggests the adaptation of arbitrary normal mutation distributions.
Choosing the mean value to be zero seems to be self-evident and is in accordancewith Beyer and Deb (2000). The current parents must be regarded as the best approximation to the optimum known so far. Using nonzero mean is equivalent to moving the
population to another place in parameter space. This can be interpreted as extrapolation. We feel extrapolation is an anti-thesis to a basic paradigm of evolutionary computation. Of course, compared to simple ESs, extrapolation should gain a significant
speed up on many test functions. But, compared to an ES with derandomized strategy
parameter control of the complete covariance matrix, the advantage from extrapolation
seems small. Correspondingly, we find algorithms with nonzero mean (Ostermeier,
1992; Ghozeil and Fogel, 1996; Hildebrand et al., 1999) unpromising.
In our opinion, any adaptation mechanism that adapts a general linear encoding/decoding has to meet the following fundamental demands:
Adaptation: The adaptation must be successful in the following sense: After an adaptation phase, progress-rates comparable to those on the (hyper-)sphere objective9
This
function
 must be realized on any convex-quadratic objective function.
must hold especially for large problem conditions (greater than , ) and nonsystematically oriented principal axes of the objective function.

  Performance: The performance (with respect to objective function evaluations)the on
(hyper-)sphere function should be comparable to the performance of the
-ES


 


, .10 A loss by a factor up towith optimal step size, where        
ten seems to be acceptable.







  

 -ES with isotropic mutation distribuInvariance: The invariance properties of thetion with respect to transformation of object parameter space and function value
should be preserved. In particular, translation, rotation, and reflection of object
parameter space (and initial point accordingly) should have no influence on strategy behavior as well as any strictly monotonically increasing, i.e., order-preserving
transformation of the objective function value.
Stationarity: To get a reasonable strategy behavior in cases of low selection pressure,some kind of stationarity condition on the strategy parameters should be fulfilled
under purely random selection. This is analogous to choosing the mean value
to be zero for the object parameter mutation. At least non-stationarity must be
quantified and judged.
From a conceptual point of view, one primary aim of strategy parameter adaptation is to become invariant (modulo initialization) against
certainof transformations
the object parameter space. For any fitness function






,
the
global
step size adaptation can yield invariance against  , . An adaptive general linear

!

!

!

!

9 That is any function  !"$#&%!'(#  , where the Hessian matrix % F H J)J is symmetric andpositive definite, ' F1HKJ and F1H . Due to numerical requirements,
the condition of % may be restricted to,

e.g., +* , . We use the expectation of J .-0/ ?+? 1( 235 4 687E ?9 1( 2:354  67E+E with ; F< as the progress measure(where  7 is the fitness value at the optimum), because it yields comparable values independent of % , ' ,measure yields values close to the common normalized progress ,  , and ; . For %D andMlarge
N G  ,Q E thismeasure =?>@  ACBEJ D F ?5G   (Beyer, 1995), where G is the distance to the optimum.
 
10 Not taking into account weighted recombination, the theoretically optimal ?  
:E -ES yields approximately *IH J* .

8Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation

encoding/decoding
can additionally yield invariance against the full rank  matrix . This invariance is a good starting point for achieving the adaptation demand.
Also, from a practical point of view, invariance is an important feature for assessing asearch algorithm (compare Section 6). It can be interpreted as a feature of robustness.
Therefore, even to achieve an additional invariance seems to be very attractive.
Keeping in mind the stated demands, we discuss two approaches to the adaptationof arbitrary normal mutation distributions in the next two sections.
3.1 A (Purely) Mutative Approach: The Rotation Angle Adaptation
Schwefel (1981) proposed a mutative approach to the adaptation of arbitrary normal
mutation distributions with zero mean, often referred to as correlated
mutations. The
distribution is parameterized by  standard deviations and   / rotation angles
that undergo a selection-recombination-mutation
scheme.
The algorithm of a
 -CORR-ES as used in Section 7.2 is briefly described.
' 
At generation  for each offspring parents are selected. The (componentwise) arithmetic mean of step sizes, angles, and object parameter vectors of the par
ents, denoted with,, and , are starting points for the mutation. Step sizes and angles read component wise











        







   








!

5

,











  5  , )   


(7)


   5  ,         /
',

 

(8)



  The random number 5 ,(7), which denotes a normal distribution  in Equation 
'
with zero mean and standard deviation
 , is only drawn once for all  .
The modulo operation ensures the angles to be in the interval /, whichis, in our experience, of minor relevance. The chosen distribution variances reflect the
typical parameter setting (see Bäck and Schwefel (1993)). The object parameter mutation reads

 






87 :9#"%$             '& )(**+






!



,'--



..

.



 .



87 :9

5



(9)

The spherical 5
is multiplied by a diagonal step size matrix deter   distribution
mined by  . The resulting
 axis-parallel mutation ellipsoid is newly oriented by successive rotations in all (i.e.,  '/ ) two-dimensional subspaces  spannedby canonical unit vectors. This complete rotation matrix is denoted with
 . This algorithm allows one to generate any  -dimensional normal distribution with zero mean
(Rudolph, 1992).
It is generally recognized that the typical values forand,, arenot sufficient for this adaptation mechanism (Bäck et al., 1997). Due to the mutative  ,approach, the parent number has presumably to scale with   . Choosing
which is roughly the number of free strategy parameters, performance on the sphereproblem declines with increasing problem dimension and becomes unacceptably poor
, . The performance demand cannot be met. This problem is intrinsic to thefor (purely) mutative approach.
The parameterization by rotation angles causes the following effects:



"







0 /

Evolutionary Computation Volume 9, Number 2

9N. Hansen and A. Ostermeier

Figure 2: Lines  of equal probability
density of two normal distributions (<) with 
 (left) and ',, (right) before and after rotation by  . The rightaxis ratios of
figure is enlarged in the  -direction; the rectangle covers the identical area in the left
and in the right. Comparing
the cross section areas of the ellipsoids, the distributions   are very similar, while the distributions with axis ratio
',,with axis ratio ofdiffer substantially.

The effective mutation strength strongly depends on the position in strategy parameter space. Mutation ellipsoids with high axis ratios are mutated much morestrongly than those with low axis ratios: Figure 2 shows two distributions withdifferent axis ratios before and after rotation by the typical variance. The resulting
cross section area, a simple measure of their similarity, differs in a wide range. For axis ratio(spheres, not shown) it is ,, . For axis ratio   (Figure 2 left)
it is about ', , while it decreases to roughly  for axis ratio
,', (Figure 2
right). The section area tends to become zero for increasing axis ratios. The mutation strength (on the rotation angles) cannot be adjusted independently of theactual position in strategy parameter space.
In principle, invariance against a new orientation of the search space is lost, because rotation planes are chosen with respect to the given coordinate system. In
simulations, we found the algorithm to be dependent even on permutations of thecoordinate axes (Hansen et al., 1995; Hansen, 2000)! Furthermore, its performance
depends highly on the orientation of the given coordinate system (compare Section
7.2). The invariance demand is not met.
Taking into account these difficulties, it is not surprising to observe that the adaptation
demand is not met. Progress rates on convex-quadratic functions with high axis ratios
can be several orders of magnitude lower than progress rates achieved on the sphereproblem (Holzheuer, 1996; Hansen, 2000, and Section 7.2).
When the typical intermediate recombination is applied to the step sizes, they increase unbounded under random selection. The systematic drift is slow, usually causes
no problems, and can even be an advantage in some situations. Therefore, we assume
the stationarity demand to be met.
Using another parameterization together with a suitable mutation operator cansolve the demands for adaptation and invariance without giving up the concept ofMSC (Ostermeier and Hansen, 1999). To satisfy the performance demand the conceptof MSC has to be modified.
3.2 A Completely Derandomized Approach: The Covariance Matrix Adaptation(CMA)
The covariance matrix adaptation (Hansen and Ostermeier, 1996) is a second-level (i.e.,completely) derandomized self-adaptation scheme. First, it directly implements the10

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation

aim of MSC to raise the probability of producing successful mutation steps: The covariance matrix of the mutation distribution is changed in order to increase the probabilityof producing the selected mutation step again. Second, the rate of change is adjusted

according to the number of strategy parameters to be adapted (by means of  in
Equation(15)). Third, under random selection, the expectation of the covariance matrix is stationary. Furthermore, the adaptation mechanism is inherently independentof the given coordinate system. In short, the CMA implements a principle component
analysis of the previously selected mutation steps to determine the new mutation distribution. We give an illustrative description of the algorithm.
Consider a special method to produce realizations of a  -dimensional normal dis%   %
tribution with zero mean.
, , span , and 5 ,   If the vectors
denotes independent ,  -normally distributed random numbers, then





5  ,    %     5  ,    % 

(10)

 !



 ,is a normally distributed random vector with zero mean.11 Choosing % , appropriate any normal distribution with zero mean can be realized.
The distribution in Equation (10) is generated by adding the “line distributions”
5  ,   % 3 5 87  % %  . If the vector % is given, the normal (line) distribution withzero mean, which produces the vector % 8withhighest probability (density) of all
7  % % the  . (The proof is simple and omitted.)normal distributions with zero mean, is 5
The covariance matrix adaptation (CMA) constructs the mutation distribution outof selected mutation steps. In place of the vectors % in Equation (10), the selected
mutation steps are used with exponentially decreasing weights. An example is shownin Figure 3, where  . The isotropic initial distribution is constructed by means of 
the unit vectors  and . At every generation  , the selected mutation step 12
%   is added to the vector tuple. All other vectors are multiplied by a factor.
According to Equation (10), after generation @  , the distribution reads



 





 



5  ,       5  ,      5  ,     %    5  ,    %  (  5  ,    %  (11)This mutation distribution tends to reproduce the mutation steps selected in thepast generations. In the end, it leads to an alignment of the distributions before and after selection, i.e., an alignment of the recent mutation distribution with the distributionof the selected mutation steps. If both distributions become alike, as under random
selection, in expectation, no further change of the distributions takes place (Hansen,
1998).
This illustrative but also formally precise description of the CMA differs in threepoints from the CMA-ES formalized in Section 5. These extensions are as follows:
Apart from the adaptation of the distribution shape, the overall variance of themutation distribution is adapted separately. We found the additional adaptationof the overall variance useful for at least two reasons. First, changes of the overallvariance and of the distribution shape should operate on different time scales. Dueto the number of parameters to be adapted, the adaptation of the covariance matrix, i.e., the adaptation of the distribution shape, must operate on a time scale of  . Adaptation of the overall variance should operate on a time scale  , because thevariance should be able to change as fast as required on simple objective functions
11 The covariance matrix of this normally distributed random vector reads 

12 

adjusts the change rate, and  L corresponds to #Evolutionary Computation Volume 9, Number 2

in Equation (15).

#  H H+H   # .
11

N. Hansen and A. Ostermeier














% 
 -,

 

 



  %  : 

 %  (%  %    % 





 
  

Figure 3: Construction of the mutation distribution in the CMA, where 0. The
initial configuration (upper left) consists of the orthogonal unit vectors  and . They 
produce an isotropic mutation distribution (dashed). The vectors %   %  ( , and %  ' 
 , while old vectors are multiplied
are added successively at generations  by -,  . The covariance matrix of the distribution after generation   reads   %   %               . The numbers             and  are eigenvalues of   ;  and are the corresponding eigenvectors with unit
length. Realizations
from the distribution can be drawn with Equation (11) or more









easily with 5 ,     5 , .

   $ 
&
 

 

 . Second, if overall variance is not adapted faster than thelike the spheredistribution shape, an (initially) small overall variance can jeopardize the search
process. The strategy “sees” a linear environment, and adaptation (erroneously)
enlarges the variance in one direction.
12

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation

The CMA is formulated for parent number 
and weighted multi-recombination, resulting in a more sophisticated computation of %   .


The technique of cumulation is applied to the vectors, that construct the distribution. Instead of the single mutation steps %   , evolution paths  are used. An
evolution path  represents a sequence of selected mutation steps. The technique of cumulation is motivated in the next section.
The adaptation scheme is formalized by means of the covariance matrix of the
distribution, because storing as many as  vectors is not practicable. That is, in every
generation  , after selection of the best search points has taken place,






1. The covariance matrix of the new distribution   is calculated due to 
andthe vector %   (with cumulation  ). In other words, the covariance matrix of the
distribution is adapted.
2. Overall variance is adapted by means of the cumulative path length control. Referring to Equation (5), %  is replaced by a “conjugate” evolution path ( in Equations (16) and (17)).
3. From the covariance matrix, the principal axes of the new distribution and theirvariances are calculated. To produce realizations from the new distribution the line distributions are added, which correspond to the principal axes of the mutation distribution ellipsoid (compare Figure 3).

Storage requirements are    . We note, that “generally, for any moderate  , this isan entirely trivial disadvantage” (Press et al., 1992). For computational and numericalrequirements, refer to Sections 5.2 and 7.1.

4

Utilizing the Evolution Path: Cumulation

The concept of MSC utilizes selection information of a single generation step. In contrast, it can be useful to utilize a whole path taken by the population over a number ofgenerations. We call such a path an evolution path of the ES. The idea of the evolution
path is similar to the idea of isolation. The performance of a strategy can be evaluatedsignificantly better after a couple of steps are taken than after a single step. It is worth
noting that both quantity and quality of the evaluation basis can improve.
Accordingly, to reproduce successful evolution paths seems more promising thanto reproduce successful single mutation steps. The former is more likely to maximize
a progress rate while the latter emphasizes the selection probability. Consequently,
we expand the idea of strategy parameter control: The adaptation should maximize
the probability to reproduce successful, i.e., actually taken, evolution paths rather thansingle mutation steps (Hansen and Ostermeier, 1996).
An evolution path contains additional selection information compared to single
steps – correlations between single steps successively taken in the generation sequence
are utilized. If successively selected mutation steps are parallel correlated (scalar product greater than zero), the evolution path will be comparatively long. If the steps are
anti-parallel correlated (scalar product less than zero), the evolution path will be comparatively short. Parallel correlation means that successive steps are going in the same
direction. Anti-parallel correlation means that the steps cancel each other out. We
assume both correlations to be inefficient. This is most obvious if the correlation/anticorrelation between successive steps is perfect. These steps can exactly be replacedwith the enlarged/reduced first step.
Evolutionary Computation Volume 9, Number 2

13

N. Hansen and A. Ostermeier




/ 

/ 

Figure 4: Two idealized evolution paths (solid) in a ridge topography (dotted). The
distributions constructed by the single steps (dashed, reduced in size) are identical.
Consequently, to maximize mutation step efficiency, it is necessary to realize longer
steps in directions where the evolution path is long –– here the same distance can be covered by fewer but longer steps. On the other side, it is appropriate to realize shorter
steps in directions where the evolution path is short. Both can be achieved using evolution paths rather than single mutation steps for the construction of the mutation distribution as described in Section 3.2.
We calculate an evolution path within an iterative process by (weighted) summation of successively selected mutation steps. The evolution path
obeys












/

 

   ./   %   





(12)


7


where ,and  2. This procedure, introduced in Ostermeier et al. (1994b),

 

is called cumulation.
/  is a normalization factor: Assuming %    and 3
in Equation (12) to be independent and identically
distributed yields 


 
independently ofand  are identical because <, . Variances of 

 
 




/  /   . If  , no cumulation takes  place, %   . The life
and 

 
span of information accumulated in  is roughly(Hansen, 1998): After about,  . That
generations
  the original information in  is reduced by the factor 
means
can roughly be interpreted as the number of summed steps.
The benefit from cumulation is shown with an idealized example. Consider thetwo different sequences of four consecutive generation steps in Figure 4. Compare any
evolution path, that is, any sum of consecutive steps in the left and in the right of the
figure. They differ significantly with respect to direction and length. Notice that thedifference is only due to the sign of vectors two and four.
Construction of a distribution from the single mutation steps (single vectors in
the figure) according to Equation (10) leads to exactly the same result in both cases(dashed). Signs
  of constructing vectors do not matter, because the vectors are multiplied by a ,  -normally distributed random number that is symmetrical about zero.
The situation changes when cumulation is applied. We focus on the left of Figure 4.
Consider a continuous sequence of the shown vectors –– i.e., an alternating sequence of' 
the two vectors and  . Then %  
and %     for , and accordingto Equation (12), the evolution path  alternately reads
















14





  /  $



/

 







/






!


/



 






/





  % &

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation

S H 
S J H S H 





+*

Figure 5: Resulting cumulation vectors  and  and distribution ellipsoidsfor analternating selection of the two vectors and  from Figure 4, left. For 
, no

 
cumulation
takes
place,
and
and(upperleft).
Shownistheresult
for
 
 



  )  ) ,  ) ,  ) ,', , respectively.
values ofat odd generation number  and






  /  $   / at even generation number. For  
















 



 



/



 



/










/





  % &

we get









/

/













 
 

/
/



and

 



After generations, the deviation from these limits is under  .
To visualize the effect of cumulation, these vectors and the relatedare distributions
shown in Figure 5 for different values of . With increasing life span
, the vectors and  become more and more similar. The distribution scales down in the horizontal
direction, increasing in the vertical direction. Already for a life span of approximately
three generations, the effect appears to be very substantial. Notice again that this effectEvolutionary Computation Volume 9, Number 2

15

N. Hansen and A. Ostermeieris due to correlations between consecutive steps. This information cannot be detectedutilizing single mutation steps separately. Because additional information is evaluated,cumulation can make adaptation of the distribution parameters faster and more reliable.
In our research, we repeatedly found cumulation to be a very useful concept. In
the following CMA-ES algorithm, we use cumulation in two places: primarily for the
cumulative path length control, where cumulation takes place in Equation (16), and thestep size is adapted from the generated evolution path  in Equation (17); second,
for the adaptation of the mutation distribution as exemplified in this section and taking
place in Equation (14).

5The-CMA-ES Algorithm

We formally describe an evolution strategy utilizing the methods illustrated in Sections
3.2 and 4, denoted by CMA-ES.13 Based on a suggestion by Ingo Rechenberg (1998),
weighted   14 recombination from the best out of individuals is applied, denoted with
 . Weighted
is a generalization of intermediate multi-recombi  recombination  , where
all recombination weights are identical. A MATLAB
nation, denoted 
implementation of the -CMA-ES algorithm is given in Appendix A.
  
The transition from generation  to  of the parameters  
,






, 
, 
, and  1
, as written in Equations (13) to (17),


7


  and   9 (unity
completely defines the algorithm. Initialization is  matrix), while the initial object parameter vector +   and the initial step size  haveto be chosen problem dependent. All vectors are assumed to be column vectors.
'(The object parameter vector   of individual @
reads














 

 







+      <  % " ! #
 
BED F






(13)





where the following apply:


 

+
   

% 







 , object parameter vector of the  individual in generation    .
      ,     , weighted mean of the best individuals of
denotes the ! best individual.
 . The index !
 generation, step size in generation  .  is the initial component wise standard
deviation of the mutation step.
 , for   '( and   ,   independent realizations of a 87 :9  %$&(')+*-,/. ) 10

2

43

2



%$





normally
distributed random vector. Components of %  are independent ,  -normally distributed.
 -matrix.  is theSymmetrical positive definite matrix of 7 (covariance9 
 .  determines
the normally distributed random vector     5 and   .        65     87     95   :7  5   :7  , which
9is a singular value decomposition of   . Initialization  .

13 The algorithm described here is identical to the algorithms in Hansen and Ostermeier (1996) setting;

, 

=<, and ><

?AC
@ B ; in Hansen and Ostermeier (1997) setting D

H+H H DFE ,  


=<

and


H+H+H  D E .
; and in Hansen (1998) setting D
14 This notation follows Beyer (1995), simplifying the ? ;  ;H CJI'E -notation for intermediate multi-recombination to ? ;H C:I'E , avoiding the misinterpretation of ; and ;KH being different numbers.

> <

16



G

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation




! 

and diagonal -diagonal matrix (step size matrix).   , for   of   are square roots of eigenvalues of the covariance matrixelements .
  -matrix (rotation matrix) that determines the coordinate sysorthogonal tem, where the scaling with   takes place. Columns of    are (defined%$as) the normalized eigenvectors of the covariance matrix  . Thediagonal element of   squared   is the corresponding eigenvalue to the%$       for  '  .  iscolumn of   ,    . That is,    .
orthogonal, i.e., The surfaces of equal probability density of   %  
are axis-parallel (hyper-)ellipsoids.   reorients these distribution
ellipsoids
tobecomecoordinate systemindependent. The covariance matrix  determines   and   apart from signsof the columns in    and permutations of columns in both matrices accordingly.
 7 (987  is      -normally
Conferring to Equation (13), notice that       5
distributed.
. For the construction of  is adaptated by means of the evolution path
the “weighted mean selected mutation step”      %   
is used. Noticethat the step size  is disregarded. The transition of   and  reads







!

!

!




















/

/


 



 



    %   B  D F   B D F   B D F 










 

  





 

,



 

$   &



(14)







(15)

where the following apply: , sum of weighted differences of points 8



$
Note that & is a symmetrical   

. Initialization





7





.

 -matrix with rank one.
 
 ,    determines
the cumulation time for , which is roughly
.
 

 



 
/
 normalizes the variance of , because  
/
 
(see
Section4).
& ') *,/. )


 



and %  & ')+*-,/. ) is chosen that under random selection   %   )have the same variance (and are identically distributed).
 



 
, with %  from Equation (13). The index
 %   & ')+*-, . )  0  2  %   3
 

%$


denotes the index of thebest individual from  
. The

2  are identical with those for + 
weights.
 
 ,   , change rate of the covariance matrix  . For    , , no change








!

!

takes place.

Equation (15) realizes the distribution change as exemplified in Section 3.2 and Figure
3.
An additional adaptation of the global step size  is necessary, taking place on aconsiderably shorter time scale. For the cumulative path length control in the CMA-ES,
Evolutionary Computation Volume 9, Number 2

17

N. Hansen and A. Ostermeiera “conjugate” evolution path  













/

 









is calculated, where scaling with 


  



 

"!



 

    %  

,





(16)

      F   BED F , BED F  B D F    BED F   B D F       / 

is omitted.





BED

,

(17)

where the following apply: , evolution path not scaled by  . Initialization   7 .











 

the cumulation time for    , which is roughly


 
 
    ,    determines/   fulfills /       (see Section 4).







.

, damping parameter, determines the possible change rate of  in thegeneration
sequence (compare Section 2, in particular Equations (3) and 7 (5)).
(95  7 (9   )  5  7
5
7 , expectation of the length of a
normally  distributed
randomvector.
A
good
approximationis  7 (Ostermeier, 1997).
)  5 /     















  



Apart from omitting the transformation with   , Equations (16) and (14) are identiforthe cumulation, instead of      %   . Undercal: In (16) we use    %   
8

7

:

9
random selection,    %   is
 -normally distributed, independently of  .
Thereby step lengths in different directions are comparable, and the expected length of
speeds up the , denoted , is well known.15 The cumulation in Equation (14) often
adaptation of the covariance matrix but is not an essential
feature(canbe set to ).
 K
For the path length control, the cumulation is essential (  must not be considerably
smaller than ). With increasing  , the lengths of single steps become more andmore similar and therefore more and more selection irrelevant.
While in Equation (13) the product matrix     only has to satisfy the equation      65     :7  , cumulative path length control in Equation (16) requires
its factorization into an orthogonal and a diagonal matrix.

 





5.1 Parameter Setting

2
,Besidessize and parent number , the strategy parameters 2    population
0,,  , and  , connected to Equations (13), (14), (15), (16), and (17), respectively,
have to be chosen.16
The default parameter settings are summarized in Table 1. In general, the selection
2
related parameters , , and 2 are comparatively uncritical and can be chosen
0
in a wide range without disturbing the adaptation procedure. We strongly recommend


2
always choosing
and the recombination weights according to 2 .
0
By definition, all weights are greater than zero. In real world applications, the defaultsettings from Table 1 are good first guesses. Only for 
the default value, does . For a quick impression, Table 2 lists the default values for a few  ,
yield






/

O
O 
replacing ?  < OOJ E J with ?  <   L  E? "E . This formulation fulfills demands analogously to those from Hansen(1998, 18f) on (17) and avoids the unattractive approximation of J . In preliminary investigations, bothvariants seem to perform equally well. If this result holds, we would prefer the latter variant.
15 Dirk Arnold (2000) suggested a remarkable simplification of Equation (17),16 On principle, the definition of parameter ; is superfluous, because ; could be implicitely defined by
setting D C HCH H C"D E* and D E CCH H H C"D  * .



18

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation



Table 1: Default parameter setting for the    
Table 2: Default





2

   5  
 



 /

 

 





    0

2






 5

/ 

7





 

 -CMA-ES.





 

!  '      

  



    values, small numbers indicate the truncated portion.

  

3

4

 



 . In the



6










 8
, +



 11
 

 15



 21







40
I  
    -CMA-ES, where    

208', 

,17 we choose
as the default. To make the strategy more robust or more explorative in caseof multimodality, can be enlarged, choosing accordingly. In particular, for unimodaland non-noisy problems,   is often most appropriate even for
  large  .
Partly for historical reasons, in this paper, we use theBased on -CMA-ES.

2 
2
our experience, for a given , an optimal choice
of
and
only
 
 achieves -CMA-ES, where0,  . The
speed-up factors less than two compared to the
optimal recombination weights depend   on the  function  to be optimized, and it remains
an open question whether the or the scheme performs better overall(using the default parameters accordingly).
If overall simulation time does not substantially exceed, 0  K  say,  generations,should be chosen smaller than in Table 1, e.g., @ ,. Increasing  ,  orand by a factor
makes the strategy more robust. If this increases the number ofneeded function evaluations, as generally to be expected, it will be typically by a factor
less than . The explorative behavior can be improved by increasing and , or up to  -CMA-ES, or  in accordance with a suffiently large initial step size.
in theThe parameter settings are discussed in more detail in the following.
where 




2



2

0






















 /



















: (Population size) In general, the equations




and

  


 


,



give a reasonable choice for . Large values linearly worsen progress rates, e.g., on
simple problems like the sphere function. Also on more sophisticated problems,
performance can linearly decrease with increasing , because the adaptation time(in generations) is more or less independent of . Values considerably smaller thanten may decline the robustness of the strategy.
  

: (Parent number) We recommend choosing
 . In the -CMA-ES, in most cases,
1993). To provide a robust,  will suffice (Beyer, 1996a; Herdy,0
strategy, large and if need be a larger ratio of
up to , are preferable (Hansenand Ostermeier, 1997). In particular for , even 
can occasionally be thebest choice.










17 The algorithm is independent from multiplication ofzero.

Evolutionary Computation Volume 9, Number 2











? D CCH H+H C"D E E with a real number greater than19

N. Hansen and A. Ostermeier : (Cumulation
 for step size) Investigations based on Hansen (1998) give strong evidence that 
can and must be chosen with respect to the problem dimension      





 



, and  . Large values,
while K for large  , most sensible values are between   i.e.,
 , slow down the possible change rate of the step size because    but can still be a good choice for certain problem instances.

 : (Damping for step size) According
to principle considerations, the damping para   K , where is near one, and    . Therefore wemeter must be  chosendefine      . Choosing  smaller than one, e.g., ,  , can yield (slightly)larger progress rates on the sphere problem. Depending on  ,  , , and , a factor

up to three can be gained. This is recommended if overall simulation time is considerably shorter   than   function evaluations. Consequently one may choose/ ,   
 . If is chosen too small, oscillating behavior of the stepsize can occur, and strategy performance may decline drastically. Larger values forlinearly slow down the possible change rate of the step size and (consequently)
make the strategy more robust.
: (Cumulation
for distribution) Basedon empirical observations, we suspect  to be a sensible choice for . Especially when long axes have to be learned,
 should be most efficient without compromising the learning of short axes.
 : (Change rate of the covariance matrix) With decreasing change rate   , reliabilityof the adaptation increases (e.g., with respect to noise) as well as adaptation time.
For   , , the CMA-ES becomes an ES with only global step size adaptation.
For strongly
problems, it can be appropriate to choose smaller change disturbed  )   with
rates like  
instead of  .

















 



5.2 Limitations and Practical Hints
Besides the limitations of any general linear encoding/decoding scheme, the limitations for the CMA, revealed so far, result from shortage of valid selection information(Hansen and Ostermeier, 1997) or numeric precision problems. The former can be due
to selection irrelevant parameters or axes, weak or distorted selection, or also, numericprecision problems. In any of these cases, a random walk on strategy parameters will
appear. In the CMA, the condition of the covariance matrix
increases unbounded
bringing search, with respect
to
the
short
principal
axes
of
,
to
a
premature end.

The initial step size  shouldbechosen
such
that
does
not tend to increasesignificantly within the initial generations. Otherwise the initially learned distribution shape can be inefficient and may have to be unlearned consuming a consider
able number of additional function evaluations. This effect of a too small  can be
avoided keeping the mutation distribution shape (i.e., the covariance matrix) initiallyconstant in case of an initially increasing . Furthermore, the prominent effect of theinitial step size on the global search performance (see Sections 7.4 and 6) should be
kept in mind.
In the practical application, a minimal
variance for the mutation steps should be"$ensured. Remember thatis thediagonal element of the step size matrix  . To
ensure a minimal variance,the standard deviation of the shortest principal axis of themutation ellipsoid, , must be restricted. If it falls below the given bound,  

20!

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation

step size should be enlarged (compare Appendix A). Problem specific knowledge
must be used for setting the bound. Numerically, a lower limit is given by the demand

 
 
  
+  
 +    ,for each unit vector . With respect to the selection
 
  
 
 
    seems desirable.
procedure, even 8    +   , 
To avoid numerical errors, a maximal condition number for  , e.g., , , should
be ensured. If the ratio between the largest and smallest eigenvalue of   exceeds
 

 
 

 
9
 
    , "/   limits the condition, , the operation  number in a reasonable way to ,  .
The numericalstability of computing eigenvectors and eigenvalues of the symmet
the built-in function eig can berical matrix usually poses no problems. In MATLAB,
used. To avoid complex solutions, the symmetry of must be explicitly enforced (compare Appendix A). In C/C++, we used the functions tred2 and tqli from Press etal. (1992), substituting float with double and setting the maximal iteration number
in tqli to ,  . Another
implementation was done in Turbo-PASCAL. For condition 
numbers up to , , we never observed numerical problems in any of these programming languages.
If the user decides to manually change object parameters during the search process – which is not recommended – the adaptation procedure (Equations (14) to (17))
must be omitted for these steps. Otherwise, the strategy parameter adaptation can beseverely disturbed.
Because change of  is comparatively slow (time scale  ), it is possible to update   and   not after every generation but after )  or even after,
 e.g.,   , generations.
This reduces the computational effort of the strategy from     to     or   , respectively. The latter corresponds to the computational effort for a fixed linear
encoding/decoding of the problem as well as for producing a realization of an (arbitrarily) normally distributed, random vector on the computer. In practical applications,
it is often most appropriate to update   and   every generation.
A simple method to handle constraints repeats the generation step, defined inEquation (13), until or at least feasible points are generated before Equations (14)to (17) are applied. If the initial mutation distribution generates sufficiently numerous feasible solutions, this method can remain sufficient during the search process dueto the symmetry of the mutation distribution. If the minimum searched for is locatedat the edge of the feasible region, the performance will usually be poor and a moresophisticated method should be used.

    

6

Test Functions

According to Whitley et al. (1996), we prefer test functions that are non-linear, nonseparable, scalable with dimension  , and resistant to simple hill-climbing. In addition,
we find it reasonable to use test functions with a comprehensible topology even if   .
This makes the interpretation of observed results possible and can lead to remarkablescientific conclusions. Table 3 gives the test functions used in this paper. The test suite18
mainly strives to test (local) convergence speed. Only functions
7 10–12 are multimodal.
Apart from functions 6–8, the global optimum
is located at .
   , where  is the object parameTo yield non-separability, we set 

ter vector according to Equation (13).
implements an angle-preserving lin
ear transformation of  , i.e., rotation and reflection of the search space. 
18 For higher dimension, even function 8  392 has a local minimum near  ?   C  C+H H+H C :E # . With thegiven initialization of

! 





and = 



 




    , the ES usually converges into the global minimum at

Evolutionary Computation Volume 9, Number 2







  .
21

N. Hansen and A. Ostermeier

Table 3: Test functions (to be minimized), where1.
        E        $     
2.
 &

3. 
              ,', ,   4.
     ,',,         




5.
      $  ,',,   & 6.
      6/    ,,      7.
   6/   ,,      



8.
      $  ,',   /     / &  9.
       
  
10.
     ,              ,
/ ,




  


  $



Function$   







)

,

B

,

$



5 











 B
5



   

)

,

,



 






,

7





,',



,',


,',







, 
,
 
,

/ ,',,/ ,',,
 
,





        "   "   /      "   "   /      "   "   /  
,7    ,  )  ,,, B where   12.            ,  )  ,,,',, B 
where  
11.







 





,









,
















 













+ 
          
77
7
  





 











           .


,





is a randomly oriented, orthonormal basis, fixed for each simulation run. Eachisrealized equally distributed
on
the
unit
(hyper-)spheresurface,dependently
drawn
so


that 
.
, if  . An algorithm to generate the basis isgiven in Figure 6. For ' 











 +  
 , we yield 2 
the initial +   the ' .  For
canonical basis, where  . For  , it holds 6  and, obviously,
a non-canonical basis, only  $remains separable.
Functions 1–5 are convex-quadratic and can be linearly transformed into the sphereproblem.19 They serve to test the adaptation demand stated in Section 3. The axis
scaling between longest and shortest principal axis of functions 3–5 is ,,', , i.e., the
problem condition is ,  . Taking into account real-world
a search strategy
 problems,should be able to handle axis ratios up to this number.
can be interpreted as aspherewith a smooth equality constraint in  direction.
 model and  %$
facilitate straight ridge topologies. The ridge points into 
direction that has to be enlarged unlimited. The sharp ridge function $is topologically invariant from the distance to the ridge peak. In particular, the local gradientis constant, independent of the distance, which is a hard feature for a local search strategy. Results can be influenced by numerical precision problems (even dependent on theimplementation ofthe algorithm). Therefore, we establish
for this  problem a minimal
step size of ,and a maximal condition number for  of ,(compare Section

 

! 

 

  
 















 







!

!

!

!

!

!

!

!

19 They can be written in the form  ? E  3 2 92 ?"! E  ?"! E#!
#! ##!  # % , where !
a full rank %$  -matrix, and the Hessian matrix %  ! # ! is symmetric and positive definite.

22

isEvolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation!   TO 

FOR

 



1. Draw components of independently ,  -normally distributed



 
 

/

(   denotes the canonical scalar product)
2.

3.





             

ROF
Figure 6: Algorithm to generate a random orthonormal basis 
et al., 1995).

  (Hansen

5.2 and Appendix A).
 , sometimes called the “banana funcThe generalized Rosenbrock functiontion,” facilitates a  bent ridge, where the global optimum, which means   is at 
 , the sum of different powers
 cannot be linearly trans . Like
6
formed into the sphere problem. Here the mis-scaling continually increases while approaching the optimum. A local search strategy presumably follows a gradually narrowing ridge.
 Functions 10–12 are multimodal.
is the generalized form of the well-known   and       are mis-scaled versions of   . The
Rastrigin function, while
 
 



.
mis-scaling between longest and shortest axis is , for'  and ,,', forFor - ', , the factor between “adjacent” axes
, respectively. We feel
 and
 is
  is a much more realistic scenario
the moderately mis-scaled Rastrigin function
than the perfectly scaled one. On the multimodal test functions, the initial step size ischosen comparatively large. With small initial step sizes, the local optimum found by
the ES is almost completely determined by the initial point that is, in each coordinate,"     equally distributed in /.
Invariance (compare also Section 3) is an important feature of a search strategy, because it facilitates the possibility of generalizing simulation results. Imagine,for   exam  / ple, translation7 invariance is not given. Then it is not sufficient to test
for, e.g., 0
. One has to test various different values for and may end up withdifferent results open to obscure interpretations. Typically, ESs are translation invariant. Simple ESs are, in addition, invariant against rotation and reflection of the search
space, i.e., they are independent of the choice of the orthonormal coordinate system  . More complex ES algorithms may be lacking this invariance (compare Sections 3.1 and 7.2).20Any evaluation of search strategies (for example, by test functions)has to take into account this important point.

   





























7



Simulation Results

Four different evolution strategies are experimentally investigated.

    -CMA-ES, where the default parametersetting from Table 1 is used apart'
from and as given below, and  for all ! . To reduce the compu



2





tational effort, the update of   and   from  for Equations (13), (14), and(16) is done every )  generations in all simulations.
20 For example, invariance against rotation is lost when discrete recombination is applied.

Evolutionary Computation Volume 9, Number 2

23

N. Hansen and A. Ostermeier I  ,',  -CORR-ES (correlated mutations), an ES with rotation angle adapta-

tion according to Section 3.1. The initialization of the rotation angles is randomly
uniform between /and , identical for all initial individuals.
 
 

the same strategy as the,   -PATH-ES,,  -CMA-ES, while  is set to9' 
zero.  
for all @ ,, and therefore,    and     can be set to9. Equations (14) and (15) become superfluous, and only cumulative path length
control takes place.











   ,  -MUT-ES, an ES with mutative adaptation of one (global) step size and in-termediate recombination on object and strategy parameters as in CORR-ES. Mu#
tation on the strategy parameter level is carried out as in Equation (1), where  is , ')    -normally distributed.
In the beginning of Section 3, we formulated four fundamental demands on analgorithm that adapt a general linear encoding in ESs. Now we pursue the question
whether the CMA-ES can satisfy these demands.
The invariance demand is met because the algorithm is formulated inherently independent of the coordinate system: All results of CMA-ES, PATH-ES, and MUT-ES
actually chosen (see Section 6), i.e., valid forare independent of the basis 
any orthonormal coordinate system. That means, these strategies are, in particular,
independent of rotation and reflection of the search  space.
 $  -CMA-ES performs on
The performance demand is also satisfied. The   
 -ES with isotropic mubetween( large) and  (6
) times slower than thetation distribution and optimal step size. The applicability of the CMA algorithm isindependent of any reasonable choice of and .
The
stationarity
demandis respected because, under
random selection, expectation of 
equals  , and expectation of   equals   (Hansen, 1998).
Simulation runs on functions 3–5 in Section 7.2 will show that the adaptation demand is met as well.
In Section 7.3, performance results on the unimodal test functions 1–9 are presented, and scaling properties are discussed. In Section 7.4, global search properties
are evaluated.









7.1 CPU-Times
We think strategy internal time consumption is of minor relevance. It seems needless to
say that results strongly depend on implementational skills and the programming language chosen. Nevertheless, to complete thewe summarize CPU-times taken
 picture,by strategy internal operations of different,  -ESs in Figure 7. The strategies werecoded in MATLAB and run uncompiled on a pentium 400 MHz processor. Differences
between MUT-ES and PATH-ES are presumably due to the more sophisticated recombination operator in MUT-ES, which allows two parent recombination independent of
the parent number . MATLAB offers relatively efficient routines for the most time consuming computations in the CMA-ES. Therefore, even for   ', , the CMA-ES needs
only , ms CPU-time per generated offspring. This time can still be reduced doing the, generations instead of every )  generations. No
update of   and   every 
efficient MATLAB routines are available for the rotation procedure in the CORR-ES. To
make CPU-times comparable, this procedure was reimplemented in C and called from
MATLAB gaining a remarkable speed up. CORR-ES is still seven to ten times slower
than CMA-ES. This may be due to a trivial implementational matter and does not makeany difference for many non-linear, non-separable search problems.





24

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation

CORR

CPU−Times [s]−110

seconds

CMA

−2

10

−3

10

MUT
PATH
−4

10

5

10

20

40

80

160

320

dimension

Figure 7: Strategy internal  CPU-times per function evaluation (i.e., per generated offspring) of different,  -ESs in seconds. The strategies were coded in MATLAB(compare text) and run on a pentium 400 MHz processor.



7.2 Testing Adaptation
In this section, single simulation
runs on the convex-quadratic functions 3–5 are shown
in comparison with runs on  $ , revealing whether the adaptation demand can be
met.
 I

,',  -CORR-ES.
We start with a discussion of the
In Figure  8, with  respect


to
the
stop
criterion,
the
best,
middle,
and
worst
out
of
runs
on
,
, and
  
  ,,  -CMA-ES are shown
are shown, where   .21 Respective results for thefor comparison. This selection scheme is not recommended for the  CMA-ES, where  (see Section 5.1, in particular Table 1 and 2). The recommended
 -CMA-ES
performs roughly ten times faster (compare also Figure 10).
First, we compare the results between the left and right column in the figure. On
the left, the axis-parallel, and therefore completely separable, versions of the functions
are used ( and   , compare Section 6). Only axis-parallel mutation ellipsoids
isare necessary to meet the adaptation demand. On the right, the basis randomly chosen anew for each run. The CORR-ES performs roughly ten to forty times
slower on the non-axis-parallel oriented versions. In accordance with previous results(Hansen et al., 1995; Holzheuer, 1996; Hansen, 2000), the CORR-ES strongly depends
on the orientation of the given coordinate system and largely exploits the separabilityof the problem. In contrast, the CMA-ES performs identically on the left and right.
Only on the axis-parallel versions of   and 
does the CORR-ES partly meet
the adaptation
demand.
Attimes,
afteranadaptation
phase,
progress rates are similar $to those on. On the non-axis-parallel
functions,theprogress
rates are worse by
  
 

a factor
between
'
,
,
(best
run
on
)
and
',

'
,
,
(worstrunon)comparedto thoseon  $ . As Rudolph (1992) pointed out, the search problem on the strategy parameters is multimodal. Even after a supposed adaptation phase, on most test functions,long phases with different progress rates can be observed. This suggests the hypothesis that the ratio between the effect of the mutation (of the angles) and the width of the



























21 With   runs, one would expect roughly 
worst run.

of any simulations to end up between the shown best and

Evolutionary Computation Volume 9, Number 2

25

N. Hansen and A. Ostermeier2

2

Cigar, Non−Axis−Parallel

Cigar, Axis−Parallel
0

0

−2

−2

−4

−4

−6

−6

−8

−8

−10
0

−10
0

1000
2000
3000
4000
generations = feval/100

210000 20000 30000 40000
generations = feval/100

2

Tablet, Axis−Parallel

Tablet, Non−Axis−Parallel

0

0

−2

−2

−4

−4

−6

−6

−8

−8

−10
0

200

400
600
800
generations = feval/100

1000

2

−10
0

2000

4000
6000
8000
generations = feval/100

10000

2

Ellipsoid, Axis−Parallel

Ellipsoid, Non−Axis−Parallel

0

0

−2

−2

−4

−4

−6

−6

−8

−8

−10
0

100

200
300
400
500
generations = feval/100

600

−10
0

2000 4000 6000 8000 10000 12000
generations = feval/100


 I

,',  -CORR-ES ( ),Figure 8:    function value vs. generation number with
where  
. The best, middle, and worst out of  runs are shown.
 K  symbols at
Bold
 ,
the lower edge indicate the mean generation number to reach. In the leftcolumn,, that is, the functions are axis-parallel oriented and completely separable. The
range of the abscissa between left and right columns is enlarged by a factor, for   (lower row) and , otherwise. For  comparison
are shown a single run on
 '
 $
( ) and respective simulations with the, ',   -CMA-ES ( ) that do not reflect
the real performance of the CMA-ES. The default
 -CMA-ES performs roughly
ten times faster. The CORR-ES does not meet the adaptation demand. Progress rates
on the non-separable
functions are worse by a factor between ,, and ,',', compared  $to those on
.







26













Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation

10
8
6
4
2
0
−2
−4
−6
−8
−10
−12
0

6
4
2
0
−2
−4
−6
−8
−10
−12
400
0

Cigar

100

200

300

generations = feval/10

6

Tablet

200

400

600

800

generations = feval/10

1000

5

4
0

2
0
−2

Ellipsoid−4−10
Ellipsoid

−6
−15

−8
−10
0

Sphere−5100 200 300 400 500 600 700
generations = feval/10

 '

0



200 400 600 800 1000 1200 1400
generations = feval/10



Figure 9: Simulation runs with the,  -CMA-ES, where   , . :     (function:   (smallest and largest
value).
 $ variance of the mutation distribution), :. The upper ten graphs (without symbol)
    (function value) from a run on
are the variances of the principal axes of the mutation distribution ellipsoid, sorted, on a logarithmic scale. Lower right: Simulation run onand   multiplied with ,  
8  $
until function value ,is reached and afterwards on  8 . The CMAES meets the adaptation demand.
Afteranadaptation
phase,
inallcases,
the progressrates are identical with those on  $
.






respective optimum differs significantly between different optima. Therefore, different
progress rates can be observed.
Finally, we stress reliability and replicability of the results with the CORR-ES. First,
the initialization of the angles is of major importance for the results. Initializing theangles with zero, which is a global optimum
for the axis-parallel functions, is almost
ten times slower on the axis-parallel   . Also, initializing all individuals with different
angles is still considerably worse than the results shown. Second,the  implemented . Using,order of applied rotations influences the results especially on 
ande.g., the reverse rotation order, results are significantly
worse.
  ,  -CMA-ES. Due to a small variance
We continue with a discussion of the
between different simulation runs (compare Figures 8 and 10), it is sufficient here
 tolook at single runs on different functions. Figure 9 shows simulation runs on 









Evolutionary Computation Volume 9, Number 2

27

N. Hansen and A. Ostermeier

 


(upper left) and(upper right), where  , , compared to a run on  $ . The
adaptation process can be clearly observed through the variances of the principle axesof the mutation distribution, where global step size is disregarded (upper ten graphs).
When the variances correspondto the scaling in the objective function, progress rates  is chosen
identical to those on  $are realized (notice again  that the basis
is due to the cumulation, which
randomly here). The shorter adaptation time on
detects the parallel correlations of single steps very efficiently.
 In the lower left of Figure 9, a simulation run on   is shown. Similar to,
but more pronounced, local adaptation occurs repeatedly on this function. Progressincreases and decreases a few times together with the step size. When the adaptationis completed, the variances are evenly spread in a range of ,  . In the lower
right of
Figure 9, the objective function is switched from   to (some multiple of)  $ , after
reaching function value ,. The mutation distribution adapts from an ellipsoid-likeshape back to the isotropic distribution. Adaptation
time and graphs of function valueand step size are very similar to those on   . In fact, from a theoretical point of view,
the algorithm must show exactly the same behavior in both adaptation cases (apartfrom stochastic
influences). Again, after the adaptation phase, progress rates identicalto those on  $are achieved in both cases.
Concluding these results, the adaptation demand is satisfied by the CMA-ES: Anyconvex-quadratic function is rescaled into the sphere function.













7.3 Testing Convergence Speed and Scaling 
In this section, we investigate the number of function
evaluations
to
reach
on
 '

test
functions
1–9,
where
.
The
-CMA-ESiscompared
to
the



',




,


,
,

 
step size adaptation takes place and which usually,  -PATH-ES, where only
 'global(slightly) outperforms the-MUT-ES. Depending on CPU-time resources, up to,, simulation runs are evaluated. Figure  $ 10 shows the results on the convex-quadratic
test functions
1–5(firstrow).
On, both strategies perform almost identically,  $while on "$   , the results are still comparable (upper left). The difference onis due to the decreasing variance of the covariance matrix supporting the adaptation of
the step size that is somewhat too slow (see above discussion of damping parameter  ,',', ,
If the axis ratio between the longest and shortest axis is in Section 
 5.1).
 , and   (upper right), CMA-ES outperforms PATH-ES by a factor
as on
,

 
, doesbetween  (  ,    , ) and almost ',,',', (  ). Only onfor 
the factor fall below ,', .
The lower row of Figure 10 shows the results on functions 6–9.  The CMA-ES works
 , the CMA-ES
quite well even on different kinds of ridge-like
topologies. Onis  – ', times faster than PATH-ES; $  on   # , the factor becomes ,,', (lower left).
 are omitted (lower right). On the latter,
Results
forthe
PATH-ESon   and
is reached after more thanfunction
evaluations ( ). On the former, the,PATH-ES does not reachbecause step size converges
tozero(as with  mutative
 
step size control). Setting the minimum
step
size
to
,
far
more
than
,
, function
 
evaluations are needed to reach.
We take an interesting look at the scaling properties with  (compare the sloping
grids in the figure). As with simple ESs in general,  PATH-ES scales linearly on mostfunctions. There are two exceptions. First, on "$   the PATH-ES scales quadratically. The axis ratio between longest and shortest principal axis of this function increases with increasing  . That is, not only problem dimension but also “problem difficulty” increases. Therefore, a simple ES scales worse than linear. This favors CMA-ES























/

 



28

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation

2

3

9

10
6

2

5

10

1
4

10

3

10

4
5
4

7

10

6

10

5

3

10

3 cigar4 tablet
5 ellipsoid4

10

1 sphere
2 schwefel

10

5

8

function evaluations

function evaluations

10

4

3

10
4

10

30

100

300

4

10

30

dimension

100

300

dimension

6

8

10

8

10

8
7

10

8

6

10

5

10

6

4

10

function evaluations

function evaluations

7

10

9
7

6

10

5

10

4

10
6 parabolic ridge
8 rosenbrock

3

10

4

10

30

100

7 sharp ridge
9 different powers

3

10

300

4

dimension

10

30

100

300

dimension

 
 '
,  -CMA-ES (solid graphs)
Figure 10:
evaluations to reachwith the Function
and the,  -PATH-ES (dashed graphs) vs. dimension  on the test functions 1–9,
where   ', ',  ', .  Shownvaluesdeviations. Sloping are mean  and standard   dotted lines indicate  8  G, ,',  ,,',  and       ,',  ,',  function evaluations. For the PATH-ES the missing result on (   , ) could not
be obtained in reasonable time (upper right) and results on $and  #  are$far above the shown area (lower right, compare text). On(upper left), both
strategies perform almost identically.



















 

scaling slightly better than PATH-ES
on $  E . Second and more surprisingly, per 
formance of the PATH-ES on
is almost independent of  . This effect is due tothe cumulative path length control that adapts here better (i.e., larger) step sizes withhigher dimension.
In general, we expected the CMA-ES to  scale quadraticallywith  . The number and   meet ourof free strategy parameters to be adapted is   .
$
expectations quite well. On, no adaptation is necessary, and CMA-ES scales
linearly with  , as to be expected of any ES.

Evolutionary Computation Volume 9, Number 2

29

N. Hansen and A. Ostermeier , even though desirable,
In contrast, the perfect linear scaling on 
andcomes as a surprise to us. On both functions one long axis has to be adapted. We
found the cumulation
to be responsible for this excellent
scaling. Without cumulation, 
the scaling on 
is similar to the scaling on. As mentionedcu 
 %$ above, the
 , andmulation
especially
supports
theadaptation
of
long
axes.
On
,
,
  "$ E 
, increasingly ordered, the scaling is between   and . Where the CMA-ES
scales worse than the PATH-ES, performances align with  because the progresssurpasses the adaptation procedure with  .
In summary, the CMA-ES substantially outperforms
the PATH-ES in dimensions

up to  , on all tested functions –– apart from  $ , where both strategies perform
almost identically. The CMA-ES always scales between  and   : Exactly linearly forthe adaptation of long axes and if no adaptation takes place, nearly quadratically if acontinuously changing topology demands persistent adaptation.







 



 



7.4 Testing Global Search Performance
In evolutionary computation, the aspect of global search is often emphasized. In contrast, we interpret ESs –– not taking into account the initial search phase –– as local searchstrategies. The population occupies a comparatively small area, and the horizon beyond this area is limited by the actual step size: Steps larger than some multiple ofthe distribution variance do virtually not appear.22 In our opinion, even the so-calledpremature convergence is often due to the lack of local convergence speed. When developing adaptation mechanisms, we mainly strove to address local convergence speed
and did not consider global search performance.
Consequently, even though the general judgment of global search performance isproblematic, it is sometimes argued that adaptation to the local topology of the objective function spoils global search properties of an algorithm. With respect to the
CMA-ES, we discuss this objection now: There is good reason, and some evidence, thateven the opposite is a more appropriate point of view: The local adaptation mechanismof the CMA-ES improves global   search properties.
,  -ESs on the generalized Rastrigin function, where  
We compare different, . Increasing population size improves the performance on this function. In contrast,
the differences between the strategies compared are not affected. Note that a smaller
initial step size worsens the performance as interpreted in Section 6. The CMA-ES iscompared to the MUT-ES, where mutative control of only global step size takes place.
PATH-ES performs similar to MUT-ES on the
Rastrigin functions.
 investigated with the CMA-ES (left) and the MUTFigure 11 shows ', simulation runs on
ES (right). Behavior of both strategies is very similar. They get trapped into local min
ima with function values between , and ,, within about ',',', and ,',, function
evaluations.
   that
Figure 12 shows ', simulation runs on the scaled Rastrigin function
should be regarded as a more realistic multimodal test problem than
. Even though
the mis-scaling between longest and shortest axis is only of a factor ten, the CMA-ES(left) and the MUT-ES (right) perform quite differently here. Function values obtainedwith MUT-ES are worse by a factor of ', than those   obtained with CMA-ES. In addition,, function evaluations compared tothese values are finally reached after about









22 In a pure mathematical sense, this is, of course, wrong. But, in contrast to any theoretical consideration,in practical applications, the finite time horizon is too short to wait for those events to occur. Even if a
distribution is chosen that facilitates large steps more often, due to a search space volume phenomenon,
these steps will virtually never produce better points if  exceeds, say, ten.

30

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation





7

10

CMA-ES on

6

10

5

function valuefunction value

5

10

4

10

3

10

4

10

3

10

10

2

2

10

10

1

0



MUT-ES on

6

10

10



7

10

1

200

400
600
generations

800

10

1000

0

200

400
600
generations

800



1000



 with the
Figure   11: ', simulation runs  on the generalized Rastrigin function,  -CMA-ES (left) and the,  -MUT-ES (right), where < , . Both strategies
perform very similar.





7

10

CMA-ES on 

MUT-ES on





6

6

10

10

function valuefunction value

5

10

4

10

4

10

3

10

2

2

10

10

,

1

10

0

200

400
600
generations

800

1000

0

200

400
600
generations

800

0 5 10
1000
4
x 10

 
  , maximal axis
Figure 12:
on the scaled Rastrigin  function
 runs   ', simulation
 
'
ratio
, , with the,  -CMA-ES (left) and the,  -MUT-ES (right), where  ', . Function values reached with MUT-ES are worse by a factor , .







  

,  function evaluations with the CMA-ES.

Figure 13 shows , simulation runs on        . The  variance of function values
reached with the CMA-ES (left) is larger than on   and     . Nevertheless, about
half of the simulation  runs end up in local optima with function values less than ,,  .
This8, function evaluations, i.e., about ten times longer than on  takes about or    . The adaptation time is now a decisive factor. With the MUT-ES, obtainedfunction values are worse by a factor ,,',, than those obtained with the CMA-ES,
and simulation roughly needs ten times the number of function evaluations (strictly
speaking, even after ,  function evaluations, final function values are not yet reached).

We found very similar results to those presented here on the Rastrigin functions in earlier investigation on a more complex, constrained multimodal test problem(EVOTECH-7, 1997).
Even though at first glance these results are surprising, there is a simple explanation why global search properties are improved. While the shape of the mutation
distribution becomes suitably adapted to the topology of the objective function, thestep size is adjusted much larger than without adaptation. A larger step size improves
Evolutionary Computation Volume 9, Number 2

31

N. Hansen and A. Ostermeier

10

10

CMA-ES on

10 

MUT-ES on

8

8

10function value

10

function value



10

6

10

6

10

4

10

4

10

2

10

2

10

, 
0

1000

2000
3000
generations

4000

5000

0

1000

2000
3000
generations

4000

0 5 10
5000
4
x 10

 
    , maximal axis
Figure 13:     ', simulation runs
 on the scaled Rastrigin function,,', , with the,  -CMA-ES (left) and the,  -MUT-ES (right), whereratio  ', . Function values reached with the MUT-ES are worse by a factor ,,',, .







the global search performance of a local search procedure. The effect of the distribution adaptation on the step size can be clearly observed in Figure 9. Variances in all
directions continually increase with the ongoing adaptation to the topology on 
between generations,', and  ', (shown are the smallest and the largest variance) and
several times on   .
Concluding this section, we observe the local adaptation in the CMA-ES to comealong with an increasing step size that is more likely to improve than to spoil global
search performance of a local search algorithm.

8Conclusion

For the application to real-world search problems, the evolution strategy (ES) is a relevant search strategy if neither derivatives of the objective function are at hand, nor differentiability and numerical accuracy can be assumed. If the search problem is expected
to be significantly non-separable, an ES with covariance matrix adaptation (CMA-ES),
as put forward in this paper, is preferable to any ES with only global or individual step
size adaptation.
We believe that there are principle limitations to the possibilities of self-adaptationin evolution strategies: To reliably adapt a significant change of the mutation distri
bution shape, at least roughly ,  function evaluations have to be done (where  isthe problem dimension). For a real-world search problem, it seems unrealistic to expect the adaptation to improve strategy behavior before, say, ',  function evaluations
(apart from global step size adaptation). A complete adaptation can take even ,,  function evaluations (compare Section 7.3). Therefore, to get always  the most out of
adaptation, CPU resources should allow roughly between ,, and ',',      function
evaluations.
One main reason for the general robustness of ESs is that the selection scheme
is solely based on a ranking of the population. The CMA-ES preserves this selection
related robustness because no additional selection information is used. Furthermore,
the CMA-ES preserves all invariance properties against transformations of the search -ES withspace and of the objective function value, which facilitates any simple
isotropic mutation distribution. Apart from initialization of object and strategy parameters, the CMA-ES yields an additional invariance against any linear transformation
32

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptationof the search space –– in contrast to any evolution strategy to our knowledge (besides
Hansen et al. (1995) and Hansen and Ostermeier (1999)). Invariance properties are ofmajor importance for the evaluation of any search strategy (compare Section 6).
The step from an ES with isotropic mutation distribution to an ES facilitating
the adaptation of arbitrary normal mutation distributions, if successfully taken, canbe compared with the step from a simple deterministic gradient strategy to a quasiNewton method. The former follows the local gradient, which in a certain sense, also
does an ES with isotropic mutation distribution. The latter approximates the inverse
Hessian matrix in an iterative sequence without acquiring additional information on
the search space. This is exactly what the CMA does for evolution strategies.
In simulations, the CMA-ES reliably approximates the inverse Hessian matrix ofdifferent objective functions. In addition, there are reported successful applications of
the CMA-ES to real-world search problems (Alvers, 1998; Holste, 1998; Meyer, 1998;
Lutz and Wagner, 1998a; Lutz and Wagner, 1998b; Olhofer et al., 2000; Bergener et al.,
2001; Cerveri et al., 2001; Igel and von Seelen, 2001; Igel et al., 2001). Consequently,
comparable to quasi-Newton methods, we expect this algorithm, or at least some quite
similar method, based on its superior performance to become state-of-the-art for theapplication of ESs to real-world search problems.

Acknowledgments
This work was supported by the Deutsche Forschungsgemeinschaft under grant
Re 215/12-1 and the Bundesministerium für Bildung und Forschung under grant
01 IB 404 A. We gratefully thank Iván Santibáñez-Koref for many helpful discussions
and persistent support of our work. In addition, we thank Christian Igel and all responding users of the CMA-ES who gave us helpful suggestions from many different
points of view.

A

CMA-ES in MATLAB% CMA-ES for non-linear function minimization% See also http://www.bionik.tu-berlin.de/user/nikofunction xmin=cmaes% Set dimension, fitness fct, stop criteria, start values...
N=10; strfitnessfct = ’cigar’;maxeval = 300*(N+2)ˆ2; stopfitness = 1e-10; % stop criteria
xmeanw = ones(N, 1); % object parameter start point (weighted mean)
sigma = 1.0; minsigma = 1e-15; % step size, minimal step size
% Parameter setting: selection,
lambda = 4 + floor(3*log(N)); mu = floor(lambda/2);
arweights = log((lambda+1)/2) - log(1:mu)’; % for recombination% parameter setting: adaptation
cc = 4/(N+4); ccov = 2/(N+2ˆ0.5)ˆ2;
cs = 4/(N+4); damp = 1/cs + 1;
% Initialize dynamic strategy parameters and constants
B = eye(N); D = eye(N); BD = B*D; C = BD*transpose(BD);pc = zeros(N,1); ps = zeros(N,1);
cw = sum(arweights)/norm(arweights);
chiN = Nˆ0.5*(1-1/(4*N)+1/(21*Nˆ2));% Generation loop
counteval = 0; arfitness(1) = 2*abs(stopfitness)+1;
while arfitness(1) > stopfitness & counteval < maxeval
Evolutionary Computation Volume 9, Number 2

33

N. Hansen and A. Ostermeier% Generate and evaluate lambda offspringfor k=1:lambda% repeat the next two lines until arx(:,k) is feasiblearz(:,k) = randn(N,1);
arx(:,k) = xmeanw + sigma * (BD * arz(:,k));% Eq.(13)arfitness(k) = feval(strfitnessfct, arx(:,k));counteval = counteval+1;
end% Sort by fitness and compute weighted mean[arfitness, arindex] = sort(arfitness); % minimization
xmeanw = arx(:,arindex(1:mu))*arweights/sum(arweights);
zmeanw = arz(:,arindex(1:mu))*arweights/sum(arweights);% Adapt covariance matrix
pc = (1-cc)*pc + (sqrt(cc*(2-cc))*cw) * (BD*zmeanw);
C = (1-ccov)*C + ccov*pc*transpose(pc);% adapt sigma
ps = (1-cs)*ps + (sqrt(cs*(2-cs))*cw) * (B*zmeanw);
sigma = sigma * exp((norm(ps)-chiN)/chiN/damp);% Eq.(14)% Eq.(15)% Eq.(16)% Eq.(17)% Update B and D from C
if mod(counteval/lambda, N/10) < 1
C=triu(C)+transpose(triu(C,1)); % enforce symmetry[B,D] = eig(C);% limit condition of C to 1e14 + 1
if max(diag(D)) > 1e14*min(diag(D))
tmp = max(diag(D))/1e14 - min(diag(D));
C = C + tmp*eye(N); D = D + tmp*eye(N);end
D = diag(sqrt(diag(D))); % D contains standard deviations now
BD = B*D; % for speed up only
end % if mod% Adjust minimal step size
if sigma*min(diag(D)) < minsigma ...
| arfitness(1) == arfitness(min(mu+1,lambda)) ...
| xmeanw == xmeanw ...
+ 0.2*sigma*BD(:,1+floor(mod(counteval/lambda,N)))sigma = 1.4*sigma;
end
end % while, end generation loopdisp([num2str(counteval) ’: ’ num2str(arfitness(1))]);
xmin = arx(:, arindex(1)); % return best point of last generation

function f=cigar(x)f = x(1)ˆ2 + 1e6*sum(x(2:end).ˆ2);

References
Alvers, M. (1998). Zur Anwendung von Optimierungsstrategien auf Potentialfeldmodelle. Berliner geowissenschaftliche Abhandlungen, Reihe B: Geophysik. Selbstverlag Fachbereich Geowissenschaften, Freie Universität Berlin, Germany.
Arnold, D. (2000). Personal communication.

34

Evolutionary Computation Volume 9, Number 2

Derandomized Self-AdaptationBäck, T. and Schwefel, H.-P. (1993). An overview of evolutionary algorithms for parameter optimization. Evolutionary Computation, 1(1):1–23.
Bäck, T., Hammel, U., and Schwefel, H.-P. (1997). Evolutionary computation: Comments on thehistory and current state. IEEE Transactions on Evolutionary Computation, 1(1):3–17.
Bergener, T., Bruckhoff, C., and Igel, C. (2001). Parameter optimization for visual obstacle detection using a derandomized evolution strategy. In Blanc-Talon, J. and Popesc, D., editors,
Imaging and Vision Systems: Theory, Assessment and Applications, Advances in Computation:
Theory and Practice. NOVA Science Books, Huntington, New York.
Beyer, H.-G. (1995). Toward a theory of evolution strategies: On the benefit of sex - thetheory. Evolutionary Computation, 3(1):81–110.

 -

Beyer, H.-G. (1996a). On the asymptotic behavior of multirecombinant evolution strategies. In
Voigt, H.-M. et al., editors, Proceedings of PPSN IV, Parallel Problem Solving from Nature, pages
122–133, Springer, Berlin, Germany.
Beyer, H.-G. (1996b). Toward a theory of evolution strategies: Self-adaptation. Evolutionary Computation, 3(3):311–347.
Beyer, H.-G. (1998). Mutate large, but inherit small! In Eiben, A. et al., editors, Proceedings of
PPSN V, Parallel Problem Solving from Nature, pages 109–118, Springer, Berlin, Germany.
Beyer, H.-G. and Deb, K. (2000). On the desired behaviors of self-adaptive evolutionary algorithms. In Schoenauer, M. et al., editors, Proceedings of PPSN VI, Parallel Problem Solving fromNature, pages 59–68, Springer, Berlin, Germany.
Cerveri, P., Pedotti, A., and Borghese, N. (2001). Enhanced evolution strategies: A novel approachto stereo-camera calibration. IEEE Transactions on Evolutionary Computation, in press.
EVOTECH-7 (1997). Evotech—Einsatz der Evolutionsstrategie in Wissenschaft und Technik,
7. Zwischenbericht. Interim report of the Fachgebiet Bionik und Evolutionstechnik der Technischen Universität Berlin under grant 01 IB 404 A of the Bundesminister für Bildung, Wissenschaft, Forschung und Technologie.
Ghozeil, A. and Fogel, D. B. (1996). A preliminary investigation into directed mutations in evolutionary algorithms. In Voigt, H.-M. et al., editors, Proceedings of PPSN IV, Parallel Problem
Solving from Nature, pages 329–335, Springer, Berlin, Germany.
Hansen, N. (1998). Verallgemeinerte individuelle Schrittweitenregelung in der Evolutionsstrategie.
Eine Untersuchung zur entstochastisierten, koordinatensystemunabhängigen Adaptation der Mutationsverteilung. Mensch und Buch Verlag, Berlin, Germany. ISBN 3-933346-29-0.
Hansen, N. (2000). Invariance, self-adaptation and correlated mutations in evolution strategies.
In Schoenauer, M. et al., editors, Proceedings of PPSN VI, Parallel Problem Solving from Nature,pages 355–364, Springer, Berlin, Germany.
Hansen, N. and Ostermeier, A. (1996). Adapting arbitrary normal mutation distributions in evolution strategies: The covariance matrix adaptation. In Proceedings of the 1996 IEEE International Conference on Evolutionary Computation, pages 312–317, IEEE Press, Piscataway, New
Jersey.
Hansen, N. and Ostermeier, A. (1997). Convergence propertiesof evolution strategies with the H , )-CMA-ES. In Zimmermann, H.derandomized covariance matrix adaptation: The (J., editor, Proceedings of EUFIT’97, Fifth European Congress on Intelligent Techniques and Soft
Computing, pages 650–654, Verlag Mainz, Aachen, Germany.
Hansen, N., Ostermeier, A., and Gawelczyk, A. (1995). On the adaptation of arbitrary normal
mutation distributions in evolution strategies: The generating set adaptation. In Eshelman,
L., editor, Proceedings of the Sixth International Conference on Genetic Algorithms, pages 57–64,
Morgan Kaufmann, San Francisco, California.
Evolutionary Computation Volume 9, Number 2

35

N. Hansen and A. Ostermeier

Herdy, M. (1993). The number of offspring as strategy parameter in hierarchically organizedevolution strategies. SIGBIO Newsletter, 13(2):2–7.
Hildebrand, L., Reusch, B., and Fathi, M. (1999). Directed mutation—a new self adaptation forevolution strategies. In Angeline, P., editor, Proceedings of the 1999 Congress on Evolutionary
Computation CEC99, pages 1550–1557, IEEE Press, Piscataway, New Jersey.
Holste, D. (1998). Modellkalibrierung am Beispiel von Kläranlagenmodellen. In Hafner, S., editor, Industrielle Anwendungen Evolutionärer Algorithmen, chapter 4, pages 37–44, Oldenbourg
Verlag, München, Germany.
Holzheuer, C. (1996). Analyse der Adaptation von Verteilungsparametern in der Evolutionsstrategie. Diploma thesis, Fachgebiet Bionik und Evolutionstechnik der Technischen Universität Berlin, Berlin, Germany.
Igel, C. and von Seelen, W. (2001). Design of a field model for early vision: A case study of
evolutionary algorithms in neuroscience. In 28th Goettingen Neurobiology Conference. In press.
Igel, C., Erlhagen, W., and Jancke, D. (2001). Optimization of neural fields models. Neurocomputing, 36(1–4):225–233.
Lutz, T. and Wagner, S. (1998a). Drag reduction and shape optimization of airship bodies. Journalof Aircraft, 35(3):345–351.
Lutz, T. and Wagner, S. (1998b). Numerical shape optimization of natural laminar flow bodies.
In Proceedings of 21st ICAS Congress, International Council of the Aeronautical Sciences andthe American Institute of Aeronautics, Paper No. ICAS-98-2,9,4.
Meyer, M. (1998). Parameteroptimierung dynamischer Systeme mit der Evolutionsstrategie.
Diploma thesis, Fachgebiet Bionik und Evolutionstechnik der Technischen Universität
Berlin, Berlin, Germany.
Olhofer, M., Arima, T., Sonoda, T., and Sendhoff, B. (2000). Optimisation of a stator blade used
in a transonic compressor cascade with evolution strategies. In Parmee, I., editor, Adaptive
Computing in Design and Manufacture (ACDM), pages 45–54. Springer Verlag, Berlin, Germany.
Ostermeier, A. (1992). An evolution strategy with momentum adaptation of the random number
distribution. In Männer, R. and Manderick, B., editors, Parallel Problem Solving from Nature,2, pages 197–206, North Holland, Amsterdam, The Netherlands.
Ostermeier, A. (1997). Schrittweitenadaptation in der Evolutionsstrategie mit einem entstochastisierten
Ansatz. Ph.D. thesis, Technische Universität Berlin, Berlin, Germany.
Ostermeier, A. and Hansen, N. (1999). An evolution strategy with coordinate system invariant
adaptation of arbitrary normal mutation distributions within the concept of mutative strategy parameter control. In Banzhaf, W. et al., editors, Proceedings of the Genetic and Evolutionary Computation Conference, GECCO-99, pages 902–909, Morgan Kaufmann, San Francisco,California.
Ostermeier, A., Gawelczyk, A., and Hansen, N. (1994a). A derandomized approach to selfadaptation of evolution strategies. Evolutionary Computation, 2(4):369–380.
Ostermeier, A., Gawelczyk, A., and Hansen, N. (1994b). Step-size adaptation based on non-local
use of selection information. In Davidor, Y. et al., editors, Proceedings of PPSN IV, Parallel
Problem Solving from Nature, pages 189–198, Springer, Berlin, Germany.
Press, W., Teukolsky, S., Vetterling, W., and Flannery, B. (1992). Numerical Recipes in C: The Art of
Scientific Computing. Second Edition. Cambridge University Press, Cambridge, England.
Rechenberg, I. (1973). Evolutionsstrategie, Optimierung technischer Systeme nach Prinzipien der biologischen Evolution. Frommann-Holzboog, Stuttgart, Germany.

36

Evolutionary Computation Volume 9, Number 2

Derandomized Self-Adaptation

Rechenberg, I. (1994). Evolutionsstrategie ’94. Frommann-Holzboog, Stuttgart, Germany.
Rechenberg, I. (1998). Personal communication.
Rudolph, G. (1992). On correlated mutations in evolution strategies. In Männer, R. and Manderick, B., editors, Parallel Problem Solving from Nature, 2, pages 105–114, North-Holland,
Amsterdam, The Netherlands.
Schwefel, H.-P. (1981). Numerical Optimization of Computer Models. Wiley, Chichester, England.
Schwefel, H.-P. (1995). Evolution and Optimum Seeking. Sixth-Generation Computer Technology
Series. John Wiley and Sons, New York, New York.
Whitley, D., Mathias, K., Rana, S., and Dzubera, J. (1996). Evaluating evolutionary algorithms.
Artificial Intelligence, 85:245–276.

Evolutionary Computation Volume 9, Number 237