Nikolaus Hansen

Anne Auger

Principled Design of Continuous Stochastic Search: From Theory to Practice

Invariance, at least according to Definition 1, can not be shown for p σ , but only for C -1 2 p c in place of p σ , or possibly for p σ with a less demanding definition of invariance, e.g. as invariance to the initial choice of the coordinate system.

Consequently, Proposition 9 must read, with corrections in magenta, Proposition 9 (Invariance under general linear transformations). If p σ k is omitted and replaced with C k -1 2 p c k for the update of step-size σ k and h σ , the CMA-ES is invariant under full rank linear transformations of the search space, that is, for each f : R n → R invariant under

The respective bijective state space transformation reads

Furthermore, for each f , the set H GL ( f ) is an equivalence class with identical algorithm trace T B (m k , σ k , C k , p c k ) for a state s and the initial state

Furthermore, covariance matrix adaptation also preserves rotation invariance (according to Definition 1) only approximately and Corollary 2 must read Corollary 2 (Rotation invariance). If p σ k is omitted and replaced with C k -1 2 p c k for the update of step-size σ k and h σ , the CMA-ES is invariant under search space rotations.

We have p σ = C -1 2 p c when C is constant over time and h σ ≡ 1. Invariance for p σ is still achieved with state space transformation (BCB T ) -1/2 BC 1/2 applied to p σ , given that C has not changed over the considered time period.

Because changes of C are generally small, replacing p σ with C -1 2 p c has usually only small effects on the performance of the algorithm.

Introduction: Top Down Versus Bottom Up

Let f : R n → R, x → f (x) be an objective or cost (or fitness) function to be minimized, where, in practice, the typical search space dimension n obeys 3 < n < 300. When properties of f are unknown a priori, an iterative search algorithm can proceed in evaluating solutions on f and so gather information for finding better solutions over time (black-box search or optimization). Good solutions have, by definition, a small f -value and evaluations of f are considered as the cost of search (remark the double entendre of the word cost for f ). The objective is, in practice, to find a Given: a cost function f , a parametrized family of distributions P(θ ), and λ ∈ N Initialize: k ← 0, set θ k Repeat while not happy Sample: x 1 , . . . , x λ ∼ P(θ k ) i.i.d. Update: θ k+1 = Update(θ k , x 1 , . . . , x λ , f (x 1 ), . . . , f (x λ )) k ← k + 1

Fig. 1 Stochastic search template good solution with the least number of function evaluations and, more rigorously, to generate a sequence x k , k = 1, 2, 3 . . . , such that f (x k ) converges fast to the essential infimum of f , denoted f * . The essential infimum f * is the largest real number such that the set of better search points {x ∈ R n : f (x) < f * } has zero volume.

In order to search in continuous spaces with even moderate dimension, some structure in the cost function needs to be exploited. For evolution strategies, the principle structure is believed to be neighborhood. Strong causality [START_REF] Rechenberg | Evolutionsstrategie, Optimierung technischer Systeme nach Prinzipien der biologischen Evolution. frommann-holzboog[END_REF]-the principle that small actuator changes have generally only small effects-and fitnessdistance correlation [START_REF] Jones | Fitness distance correlation as a measure of problem difficulty for genetic algorithms[END_REF]-a statistical perspective of the same concept-are two ways to describe the structure that evolution strategies are based upon. In contrast to Chapters ??, ?? and ?? of this volume [?], in this chapter we do not introduce an a priori assumption on the problem class we want to address, that is, we do not assume any structure in the cost function a priori. However, we use two ideas that might imply the exploitation of neighborhood: we assume that the variances of the sample distribution exist, and we encourage consecutive iteration steps to become, under a variable metric, orthogonal (via step-size control). Empirically, the latter rather reduces the locality of the algorithm: the step-sizes that achieve orthogonality are usually large in their stationary condition. We conjecture therefore that either the mere existence of variances and/or the "any-time" approach that aims to improve in each iteration, rather than only in a final step, implies already the exploitation of a neighborhood structure in our context.

In order to solve the above introduced search problem on f , we take a principled stochastic (or randomized) approach. We first sample points from a distribution over the search space with density p(.|θ ), we evaluate the points on f and finally update the parameters θ of the distribution. This is done iteratively and defines a search procedure on θ as depicted in Fig. 1. Indeed, the update of θ remains the one and only crucial element-besides the choice of p (and λ ) in the first place. Consequently, this chapter is entirely devoted to the question of how to update θ .

Before we proceed, we note that under some mild assumptions on p, and for any increasing transformation g : R → R (in particular also for the identity), the minimum of the function

θ → E(g( f (x))|θ ) (1)
coincides with the minimum of f (the expectation E is taken under the sample distribution p, given parameters θ ). The optimal distribution is entirely concentrated in the arg min of f . In black-box search, we do not want (and are not able) to impose strong regularity conditions on the unknown function f . However, we have entire control over p. This seems an excellent justification for a randomized approach to the original black-box search problem. We sketch two approaches to solve (1) 1 .

The Top Down Way

We might chose p being "sufficiently smooth" and conduct a gradient descent,

θ k+1 = θ k -η∇ θ E( f (x)|θ ) with η > 0 . (2)
We are facing two problems with Equation ( 2). On the one hand, we need to compute ∇ θ E( f (x)|θ ). On the other hand, the gradient ∇ θ strongly depends on the specifically chosen parameterization in θ . The unique solution to the second problem is the natural gradient. The idea to use the natural gradient in evolution strategies was coined in [START_REF] Wierstra | Natural evolution strategies[END_REF] and elegantly pursued in [START_REF] Glasmachers | Exponential natural evolution strategies[END_REF]. The natural gradient is unique, invariant under reparametrization and in accordance with the KL-divergence or relative entropy, the informational difference measure between distributions. We can reformulate (2) using the natural gradient, denoted ∇, in a unique way as

θ k+1 = θ k -η ∇E( f (x)|θ ) . (3)
We can express the natural gradient in terms of the vanilla gradient ∇ θ , using the Fisher information matrix, as ∇ = F -1 θ ∇ θ . Using the log-likelihood trick, ∇ θ p = (p/p)∇ θ p = p∇ θ log p we can finally, under mild assumption on p, re-arrange (3) into

θ k+1 = θ k -ηE( f (x) expensive "controlled" F -1 θ ∇ θ log p(x|θ )) . (4)
In practice, the expectation in (4) can be approximated/replaced by taking the average over a (potentially small) number of samples, x i , where computing f (x i ) is assumed to be the costly part. We will also choose p such that we can conveniently sample from the distribution and that the computation (or approximation) of F -1 θ ∇ θ log p is feasible. The top down way of (3) and ( 4) is an amazingly clean and principled approach to stochastic black-box optimization.

The Bottom Up Way

In this chapter, we choose a rather orthogonal approach to derive a principled stochastic search algorithm in the R n . We take a scrutinizing step-by-step road to construct the algorithm based on a few fundamental principles-namely maximal entropy, unbiasedness, maintaining invariance, and, under these constraints, exploiting all available information and solving simple functions reasonably fast.

Surprisingly, the resulting algorithm arrives at (3) and (4): Equations ( 12) and (51) implement (3) in the manifold of multivariate normal distributions under some monotone transformation of f [START_REF] Akimoto | Bidirectional relation between CMA evolution strategies and natural evolution strategies[END_REF][START_REF] Arnold | Information-geometric optimization algorithms: A unifying picture via invariance principles[END_REF]

(let η = 1, c 1 = c ε = 0, c µ = σ k = 1)
. The monotone transformation is driven by an invariance principle. In both ways, top down and the bottom up, the same, well-recognized stochastic search algorithm covariance matrix adaptation evolution strategy (CMA-ES) emerges. Our scrutinizing approach however reveals additional aspects that are consistently useful in practice: cumulation via an evolution path, step-size control, and different learning rates η for different parts of θ . These aspect are either well hidden by (4)2 or can hardly be derived at all (cumulation). On the downside, the bottom up way is clearly less appealing.

The following sections will introduce and motivate the CMA-ES step-by-step. The CMA-ES samples new solutions from a multivariate normal distribution and updates the parameters of the distribution, namely the mean (incumbent solution), the covariance matrix and additionally a step-size in each iteration, utilizing the franking of the sampled solutions. We formalize the different notions of invariance as well as the maximum likelihood and stationarity properties of the algorithm. A condensed final transcription of the algorithm is provided in the appendix. For a discussion under different perspectives, the reader is referred to [START_REF] Hansen | The CMA evolution strategy: A tutorial[END_REF][START_REF] Hansen | The CMA evolution strategy: a comparing review[END_REF][START_REF] Hansen | Completely derandomized self-adaptation in evolution strategies[END_REF].

Sampling with Maximum Entropy

We start by sampling λ (new) candidate solutions x i ∈ R n , obeying a multivariate normal (search) distribution

x i ∼ m k + σ k × N i (0, C k ) for i = 1, . . . , λ , (5)
where k = 0, 1, 2, . . . is the time or iteration index and m k ∈ R n , σ k > 0, and N (0, C) denotes a multivariate normal distribution with zero mean and covariance matrix C, ∼ denotes equality in distribution. For convenience, we will sometimes omit the iteration index k. New solutions obey a multivariate normal distribution with expectation m and covariance matrix σ 2 × C. Sets of equal density-that is, lines or surfaces in 2 or 3-D respectively-are ellipsoids centered about the mean and modal value m. Figure 2 shows 150 sampled points from a standard (2-variate) normal distribution, N (0, I). The circle depicts the one-σ equal density line, the center of the circle is the mean and modal value at zero. In general, lines of equal density (level sets) are ellipsoids. The probability to sample a point outside the dashed box is close to 1

-(1 -2 × 0.0015) 2 ≈ 1/170 x 1
Given mean, variances and covariances of a distribution, the chosen multivariate normal distribution has maximum entropy and-without any further knowledgesuggests itself for randomized search. We explain (5) in more detail.

• The distribution mean value, m, is the incumbent solution of the algorithm: it is the current estimate for the global optimum provided by the search procedure.

The distribution is point symmetrical about the incumbent. The incumbent m is (usually) not evaluated on f . However, it should be evaluated as final solution in the last iteration. • New solutions are obtained by disturbing m with the mutation distribution

N 0, σ 2 C ≡ σ × N (0, C) , (6)
where the equivalence holds by definition of N (., .). The parameter σ > 0 is a step-size or scale parameter and exists for notational convenience only. The covariance matrix C has n 2 +n 2 degrees of freedom and represents a full quadratic model. The covariance matrix determines the shape of the distribution, where level-sets of the density are hyper-ellipsoids (confer to [START_REF] Hansen | The CMA evolution strategy: A tutorial[END_REF][START_REF] Hansen | The CMA evolution strategy: a comparing review[END_REF] for more details). On convex quadratic cost functions, C will closely align with the inverse Hessian of the cost function f (up to a scalar factor). The matrix C defines a variable neighborhood metric. The above said suggests that using the maximum entropy distribution with finite variances implies the notion, and underlines the importance of neighborhood.

The initial incumbent m 0 needs to be provided by the user. The algorithm has no preference for any specific value and its operations are invariant to the value of m 0 (see translation invariance in Section 4).

Equation ( 5) implements the principle of stationarity or unbiasedness, because the expected value of ( 6) is zero. Improvements are not a priori made by construction, but only after sampling by selection. In this way, the least additional assumptions are built into the search procedure.

The number of candidate solutions sampled in (5) cannot be entirely derived from first principles. For small λ ≫ n the search process will be comparatively local and the algorithm can convergence fast. Only if previously sampled search points are considered, λ could be chosen to its minimal value of one-in particular if the best so-far evaluated candidate solution is always retained. We tend to disregard previous samples entirely (see below). In this case, a selection must take place between λ ≥ 2 new candidate solutions. Because the mutation distribution is unbiased, newly sampled solutions tend to be worse than the previous best solution and in practice λ ≥ 5 is advisable. 3 On the other hand, for large λ ≫ n, the search becomes more global and the probability to approach the desired, global optimum on multimodal functions is usually larger. On the downside, more function evaluations are necessary to closely approach an optimum even on simple functions.

Consequently, a comparatively successful overall strategy runs the algorithm first with a small population size, e.g. the default λ = 4 + ⌊3 ln n⌋, and afterwards conducts independent restarts with increasing population sizes (IPOP) [START_REF] Auger | A restart CMA evolution strategy with increasing population size[END_REF].

After we have established the sampling procedure using a parameterized distribution, we need to determine the distribution parameters which are essential to conduct efficient search. All parameters depend explicitly or implicitly on the past and therefore are described in their update equations.

Exploiting the Objective Function

The pairs (x i , f (x i )) i=1,...,λ , provide the information for choosing a new and better incumbent solution m k+1 as well as the new distribution covariance matrix σ 2 C. Two principles are applied.

1. Old information is disregarded. There are a few reasons to believe that old information can or should be disregarded.

• The given (n 2 + 3n)/2 distribution parameters, m and σ 2 × C, should already capture all necessary previous information. Two additional state variables, the search paths p σ , p c ∈ R n , will provide another 2n parameters. Theoretical results suggests that only slight improvements can be made by storing and using (all) previously sampled candidate solutions [START_REF] Teytaud | Lower bounds for evolution strategies using VC-dimension[END_REF][START_REF] Teytaud | General lower bounds for evolutionary algorithms[END_REF], given rank-based selection. 3 In the (µ, λ )-ES, only the µ best samples are selected for the next iteration. Given µ = 1, a very general optimality condition for λ states that the currently second best solution must resemble the f -value of the previous best solution [START_REF] Hansen | Sizing the population with respect to the local progress in (1, λ )-evolution strategies-a theoretical analysis[END_REF]. Consequently, on any linear function, λ = 2 and λ = 3 are optimal [START_REF] Hansen | Sizing the population with respect to the local progress in (1, λ )-evolution strategies-a theoretical analysis[END_REF][START_REF] Schwefel | Numerical Optimization of Computer Models[END_REF]. On the sphere function [START_REF] Hansen | Sizing the population with respect to the local progress in (1, λ )-evolution strategies-a theoretical analysis[END_REF], λ = 5 is optimal [START_REF] Rechenberg | Evolutionsstrategie, Optimierung technischer Systeme nach Prinzipien der biologischen Evolution. frommann-holzboog[END_REF]. On the latter, also λ ≈ 3.7µ can be shown optimal for µ ≥ 2 and equal recombination weights [START_REF] Beyer | The Theory of Evolution Strategies[END_REF], compare [START_REF] Hansen | The CMA evolution strategy: A tutorial[END_REF]. For λ < 5, the original strategy parameter setting for CMA-ES has been rectified in [START_REF] Brockhoff | Mirrored sampling and sequential selection for evolution strategies[END_REF], but only mirrored sampling leads to satisfactory performance in this case [START_REF] Brockhoff | Mirrored sampling and sequential selection for evolution strategies[END_REF].

• Convergence renders previously sampled solutions rather meaningless, because they are too far away from the currently focused region of interest. • Disregarding old solutions helps to avoid getting trapped in local optima.

• An elitist approach can be destructive in the presence of noise, because a supersolution can stall any further updates. Under uncertainties, any information must be used with great caution.

2. Only the ranking of the better half of the new candidate solutions is exploited. Function values are discarded as well as the ranking of the worse half of the newly sampled points. Specifically, the function f enters the algorithm only via the indices i : λ for i = 1, . . . , µ, in that (serving as definition for i : λ )

f (x 1:λ ) ≤ f (x 2:λ ) ≤ • • • ≤ f (x λ :λ ) (7)
is satisfied. We choose µ = ⌊λ /2⌋, because a. on a linear function in expectation the better half of the new solutions improve over m k and for the same reason b. on the quadratic sphere function only the better half of the new solutions can improve the performance, using positive recombination weights (see [START_REF] Hansen | The CMA evolution strategy: A tutorial[END_REF] below). For the remaining solutions, x i :λm k needs to enter with a negative prefactor [START_REF] Arnold | Weighted multirecombination evolution strategies[END_REF].

We feel that using worse points to make predictions for the location of better points might make a too strong assumption on the regularity of f in general. Indeed, optimization would be a much easier task if outstandingly bad points would allow generally valid implications on the location of good points, because bad points are generally easy to obtain.

On the highly symmetrical, isotropic sphere model, using the worse half points with the same importance than the better half points for calculating the new incumbent can render the convergence two times faster [START_REF] Arnold | Optimal weighted recombination[END_REF][START_REF] Arnold | Weighted multirecombination evolution strategies[END_REF]. In experiments with CMA-ES, we find the factor to be somewhat smaller and obtain very similar results also on the isotropic, highly multimodal Rastrigin function. On most anisotropic functions we observe performance degradations and also failures in rare cases and with noise. The picture though is more encouraging for a covariance matrix update with negative samples, as discussed below.

Because only the f -ranked solution points (rather than the f -values) are used, we denote the f -ranking also as (rank-based) selection. The exploitation of available information is quite conservative, reducing the possible ways of deception. As an additional advantage, function values do not need to be available (for example, when optimizing a game playing algorithm, a passably accurate selection and ranking of the µ best current players suffices to proceed to the next iteration). This leads to a strong robustness property of the algorithm: invariance to order-preserving transformations, see next section. The downside of using only the f -ranking is, that the possible convergence speed cannot be faster than linear [START_REF] Auger | Reconsidering the progress rate theory for evolution strategies in finite dimensions[END_REF][START_REF] Jägersküpper | Lower bounds for randomized direct search with isotropic sampling[END_REF][START_REF] Teytaud | Lower bounds for evolution strategies using VC-dimension[END_REF].

s k ✲ A f s k+1 T f →h (s k ) ✲ A h T f →h (s k+1 ) T f →h T -1 f →h ❄ ✻ T f →h ❄ Fig.

Invariance

We begin with a general definition of invariance of a search algorithm A . In short, invariance means that A does not change its behavior under exchange of f with an equivalent function h ∈ H ( f ), in general conditionally to change of the initial conditions.

Definition 1 (Invariance). Let H a be mapping from the set of all functions into its power set, H :

{R n → R} → 2 {R n →R} , f → H ( f ).
Let S be the state space of the search algorithm, s ∈ S and A f : S → S an iteration step of the algorithm under objective function f . The algorithm A is invariant under H (in other words: invariant under the exchange of f with elements of H ( f )) if for all f ∈ {R n → R}, there exists for all h ∈ H ( f ) a bijective state space transformation T f →h : S → S such that for all states s ∈ S

A h • T f →h (s) = T f →h • A f (s) , (8)
or equivalently

A h (s) = T f →h • A f • T -1 f →h (s) . (9)
If T f →h is the identity for all h ∈ H ( f ), the algorithm is unconditionally invariant under H . For randomized algorithms, the equalities hold almost surely, given appropriately coupled random number realizations, otherwise in distribution. The set of functions H ( f ) is an invariance set of f for algorithm A .

The simplest example where unconditional invariance trivially holds is H : f → { f }. Any algorithm is unconditionally invariant under the "exchange" of f with f . The idea of invariance is depicted in the commutative diagram in Fig. 3. The two possible paths from the upper left to the lower right are reflected in Equation [START_REF] Auger | Experimental comparisons of derivative free optimization algorithms[END_REF].

Equation ( 9) implies (trivially) for all k ∈ N that

A k h (s) = T f →h • A k f • T -1 f →h (s) , (10)
where A k (s) denotes k iteration steps of the algorithm starting from s. Equation [START_REF] Brockhoff | Mirrored sampling and sequential selection for evolution strategies[END_REF] reveals that for all h ∈ H ( f ), the algorithm A optimizes the function h with initial state s just like the function f with initial state T -1 f →h (s). In the lucky scenario, T f →h is the identity and A behaves identical on f and h. Otherwise, first s must be moved to T -1 f →h (s), such that after an adaptation phase any function h is optimized just like the function f . This is particularly attractive, if f is the easiest function in the invariance class. The adaptation time naturally depends on the distance between s and T -1 f →h (s). We give the first example of unconditional invariance to order-preserving transformations of f . Proposition 1 (Invariance to order-preserving transformations). For all strictly increasing functions g : R → R and for all f : R n → R, the CMA-ES behaves identical on the objective function x → f (x) and the objective function x → g( f (x)). In other words, CMA-ES is unconditionally invariant under

H monoton : f → {g • f | g is stricly increasing} . ( 11
)
Additionally, for each f : R n → R, the set of functions H monoton ( f )-the orbit of f -is an equivalence class of functions with indistinguishable search trace.

Proof idea. Only the f -ranking of solutions is used in CMA-ES and g does not change this ranking. We define the equivalence relation as f ∼ h iff ∃g strictly increasing such that f = g • h. Then, reflexivity, symmetry and transitivity for the equivalence relation ∼ can be shown elementarily recognizing that the identity and g -1 and compositions of strictly increasing functions are strictly increasing.

⊓ ⊔

The CMA-ES depends only on the sub-level sets {x | f (x) ≤ α} for α ∈ R. The monotonous transformation g does not change the sub-level sets, that is

{x | g( f (x)) ≤ g(α)} = {x | f (x) ≤ α}.

Update of the Incumbent

Given the restricted usage of information from the evaluations of f , the incumbent is generally updated with a weighted mean of mutation steps

m k+1 = m k + c m µ ∑ i=1 w i (x i:λ -m k ) (12) with µ ∑ i=1 |w i | = 1, w 1 ≥ w 2 • • • ≥ w µ , 0 < c m ≤ 1 . ( 13
)
The question of how to choose optimal weight values w i is pursued in [START_REF] Arnold | Weighted multirecombination evolution strategies[END_REF] and the default values in Table 2 of the appendix approximate the optimal positive values on the infinite dimensional sphere model. As discussed above, we add the constraints

w µ > 0 and µ ≤ λ /2 , (14)
while the formulation with ( 12) also covers more general settings. Usually, we set the learning rate c m = 1 and the computation of the new incumbent simplifies to

m k+1 = µ ∑ i=1 w i x i:λ . (15)
A learning rate of one seems to be the largest sensible setting. A value larger than one should only be advantageous, if σ k is too small and implies that the step-size heuristic should be improved. Very small σ k together with c m ≫ 1 resemble a classical gradient descent scenario.

The amount of utilized information can be quantified via the variance effective selection mass, or effective µ

µ eff = µ ∑ i=1 w 2 i -1 , (16)
where we can easily derive the tight bounds 1 < µ eff ≤ µ. Usually, a weight setting with µ eff ≈ λ /4 is appropriate. Given µ eff , the specific choice of the weights is comparatively uncritical. The presented way to update the incumbent using a Weighted mean of all µ selected points gives raise for the name (µ/µ W , λ )-CMA-ES.

Proposition 2 (Random ranking and stationarity of the incumbent). Under (pure) random ranking, m k follows an unbiased random walk

m k+1 ∼ m k + σ k √ µ eff N (0, C k ) (17)
and consequently

E(m k+1 |m k ) = m k . ( 18
)
Pure random ranking means that the index values i : λ ∈ {1, . . . , λ } do not depend on x 1 , . . . , x λ , for all i = 1, . . . , λ , for example when f (x) is a random variable with a density and does not depend on x, or when i : λ is set to i.

Proof idea. Equation ( 17) follows from Equations ( 5), ( 12) and ( 16), and (18) follows because EN (0, C) = 0 by definition.

⊓ ⊔

The proposition affirms, that only selection ( f -ranking) can induce a biased movement of the incumbent m.

Proposition 3 (Maximum likelihood estimate of the mean). Given x 1:λ , . . . , x µ:λ , the incumbent m k+1 maximizes, independent of the positive definite matrix C, the weighted likelihood

m k+1 = arg max m∈R n µ ∏ i=1 p w i N (x i:λ | m) , (19)
where p w

i N (x | m) = (p N (x | m)) w i and p N (x | m)
denotes the density of N (m, C) at point x, or equivalently the weighted log-likelihood

m k+1 = arg max m∈R n µ ∑ i=1 w i × log p N (x i:λ | m) , (20)
Proof idea. We exploit the one-dimensional normal density and the fact that the multivariate normal distribution, after a coordinate system rotation, can be decomposed into n independent mariginal distributions.

⊓ ⊔

Finally. we find translation invariance, a property that every continuous search algorithm should enjoy.

Proposition 4 (Translation invariance). The CMA-ES is translation invariant, that is, invariant under

H trans : f → {h a : x → f (x -a) | a ∈ R n } , (21)
with the bijective state transformation, T f →h a , that maps m to m + a (compare Figure 3). In other words, the trace of m k + a is the same for all functions h a ∈ H trans .

Proof idea. We consider Fig. 3: an iteration step with state (m k , σ k , C k , . . . ) using cost function x → f (x) in the upper path is equivalent with an iteration step with state (m k + a, σ k , C k , . . . ) using cost function h a : x → f (xa) in the lower path.

⊓ ⊔

Translation invariance, meaning also that m km 0 does not depend on m 0 , is a rather indispensable property for a search algorithm. Nevertheless, because m k depends on m 0 , a reasonable proposition for m 0 , depending on f , is advisable.

Step-Size Control

Step-size control aims to make a search algorithm adaptive to the overall scale of search.

Step-size control allows for fast convergence to an optimum and serves to satisfy the following basic demands on a search algorithm.

1. Solving linear functions, like f (x) = x 1 . On linear functions we desire a geometrical increase of the f -gain f (m k )f (m k+1 ) with increasing k. 2. Solving the simplest convex-quadratic function, the sphere function

f (x) = n ∑ i=1 (x i -x * i ) 2 = x -x * 2 , (22)
fast. We desire

m k -x * m 0 -x * ≈ exp -c k n , (23)
such that c ≪ 0.02 min(n, λ ), because c ≈ 0.25 λ is the optimal value which can be achieved with optimal step-size and optimal positive weights for λ ≫ n (c ≈ 0.5 λ can be achieved using also negative weights for x i :λm k in [START_REF] Hansen | The CMA evolution strategy: A tutorial[END_REF], see [START_REF] Arnold | Weighted multirecombination evolution strategies[END_REF]). The optimal step-size changes when approaching the optimum.

Additionally, step-size control will provide scale invariance, as explicated below. Unfortunately, step-size control can hardly be derived from first principles and therefore relies on some internal model or some heuristics. Line-search is one such heuristic that decides on the realized step length after the direction of the step is given. Surprisingly, a line-search can gain very little over a fixed (optimal) step length given in each iteration [START_REF] Jägersküpper | Lower bounds for hit-and-run direct search[END_REF]. Recent theoretical results even seem to indicate that in the limit for n → ∞ the optimal progress rate cannot be improved at all by a cost-free ray search on a half-line (given positive weights) or by a line search otherwise [START_REF] Jebalia | [END_REF]. A few further heuristics for step-size control are well-recognized.

1. Controlling the success rate of new candidate solutions, compared to the best solution seen so far (1/5-th success rule) [START_REF] Rechenberg | Evolutionsstrategie, Optimierung technischer Systeme nach Prinzipien der biologischen Evolution. frommann-holzboog[END_REF][START_REF] Schumer | Adaptive step size random search[END_REF]. 2. Sampling different candidate solutions with different step-sizes (self-adaptation) [START_REF] Rechenberg | Evolutionsstrategie, Optimierung technischer Systeme nach Prinzipien der biologischen Evolution. frommann-holzboog[END_REF][START_REF] Schwefel | Numerical Optimization of Computer Models[END_REF]. Selected solutions also retain their step-size. 3. Testing different step-sizes by conducting additional test steps in direction m k+1m k , resembling a rudimentary line-search (two point adaptation) [START_REF] Hansen | CMA-ES with two-point step-size adaptation[END_REF][START_REF] Salomon | Accelerating backpropagation through dynamic selfadaptation[END_REF]. 4. Controlling the length of the search path, taken over a number of iterations (cumulative step-size adaptation, CSA, or path length control) [START_REF] Ostermeier | Step-size adaptation based on non-local use of selection information[END_REF].

In our context, the last two approaches find reasonable values for σ in simple test cases (like ridge topologies). We use cumulative step-size adaptation here. The underlying design principle is to achieve perpendicularity of successive steps. Perpendicularity is measured using an evolution path and a variable metric.

Conceptually, an evolution path, or search path, of length j is the vector

m k -m k-j , (24)
that is, the total displacement of the mean during j iterations. For technical convenience, and in order to satisfy the stationary condition [START_REF] Hansen | Benchmarking a weighted negative covariance matrix update on the bbob-2010 noiseless testbed[END_REF], we compute the search path, p σ , in an iterative momentum equation with the initial path p σ 0 = 0 as

p σ k+1 = (1 -c σ ) p σ k + c σ (2 -c σ )µ eff C k -1 2 m k+1 -m k σ k . ( 25
)
The factor 1c σ > 0 is the decay weight and 1/c σ ≈ n/3 is the backward time horizon-after 1/c σ iterations about 1exp(-1) ≈ 63% of the information has been replaced; C k

-1
2 is the positive symmetric square root 4 of C k -1 . The remaining 4 The positive symmetric square root satisfies C k

-1 2 C k -1 2 = C k -1
, has only positive eigenvalues and is unique.

Fig. 4 Schematic depiction of three evolution paths in the search space (each with six successive steps of m k ). Left: single steps cancel each other out and the evolution path is short. Middle: steps are "on average orthogonal". Right: steps are positively correlated and the evolution path is long. The length of the path is a good indicator for optimality of the step-size factors are, without further degree of freedom, chosen to guaranty the stationarity

p σ k ∼ N (0, I) for k = 1, 2, 3, . . . , (26)
given p σ 0 ∼ N (0, I) and pure random ranking of x i:λ in all preceding time steps. The length of the evolution path is used to update the step-size σ either following [START_REF] Jastrebski | Improving evolution strategies through active covariance matrix adaptation[END_REF]

σ k+1 = σ k × exp c σ d σ p σ k+1 2 -n 2n (27)
or via

σ k+1 = σ k × exp c σ d σ p σ k+1 E N (0, I) -1 , (28)
where d σ ≈ 1. The step-size increases/decreases iff p σ k+1 2 or p σ k+1 is larger/smaller than its expected value. Equation ( 27) is more appealing and easier to analyze, but [START_REF] Jägersküpper | Lower bounds for randomized direct search with isotropic sampling[END_REF] might have an advantage in practice. In practice, also an upper bound to the argument of exp is sometimes useful.

Figure 4 depicts the idea of the step-size control schematically.

• If steps are positively correlated, the evolution path tends to be long (right picture). A similar trajectory could be covered by fewer but longer steps and the step-size is increased. • If steps are negatively correlated they tend to cancel each other out and the evolution path is short (left picture). Shorter steps seem more appropriate and the step-size is decreased. • If the f -ranking does not affect the length of the evolution path, the step-size is unbiased (middle picture).

We note two major postulates related to step-size control and two major design principles of the step-size update.

Postulate 1 (Conjugate steps) Successive iteration steps should be approximately C -1 -conjugate, that is, orthogonal with respect to the inner product (and metric) defined by C -1 .

As a consequence of this postulate, we have used perpendicularity as optimality criterion for step-size control.

If steps are uncorrelated, like under random selection, they indeed become approximately [START_REF] Hansen | The CMA evolution strategy: a comparing review[END_REF]. This means the steps are orthogonal with respect to the inner product defined by C -1 and therefore orthogonal in the coordinate system defined by C. In this coordinate system, the coordinate axes, where the independent sampling takes place, are eigenvectors of C. Seemingly uncorrelated steps are the desired case and achieved by using [START_REF] Hansen | Completely derandomized self-adaptation in evolution strategies[END_REF].

C -1 -conjugate, that is (m k+1 -m k ) T C -1 m k -m k-1 ≈ 0, see
C -1/2 in
In order to better understand the following assertions, we rewrite the step-size update [START_REF] Jägersküpper | Lower bounds for randomized direct search with isotropic sampling[END_REF], only using an additive update term,

log σ k+1 = log σ k + c σ d σ p σ k+1 E N (0, I) -1 . (29)
First, in accordance with our stationary design principle, we establish a stationarity condition on the step-size.

Proposition 5 (Stationarity of step-size). Given pure random ranking and p σ 0 ∼ N (0, I), the quantity log σ k performs an unbiased random walk (see Equation [START_REF] Jastrebski | Improving evolution strategies through active covariance matrix adaptation[END_REF]). Consequently, the step-size obeys the stationarity condition

E(log σ k+1 |σ k ) = log σ k . (30)
Proof idea. We analyze the update equations ( 29) and [START_REF] Hansen | Completely derandomized self-adaptation in evolution strategies[END_REF].

⊓ ⊔

Postulate 2 (Behavior on linear functions [START_REF] Hansen | An analysis of mutative σ -self-adaptation on linear fitness functions[END_REF]) On a linear function, the dispersion of new candidate solutions should increase geometrically fast in the iteration sequence, that is, linearly on the log scale. Given σ β k as dispersion measure with β > 0, we can set w.l.o.g. β = 1 and demand for some α > 0

E(log σ k+1 |σ k ) ≥ log σ k + α . (31)
The CMA-ES satisfies the postulate for some k 0 and all k ≥ k 0 , because on a linear function the expected length of the evolution path increases monotonically. We reckon that k 0 ∝ 1/c σ . Finally, we investigate the more abstract conception of scale invariance as depicted in Fig. 5.

Proposition 6 (Scale invariance). The CMA-ES is invariant under

H scale : f → {h α : x → f (x/α) | α > 0} ( 32
)
with the associated bijective state space transformation That means for all states

(m k , σ k , C k , . . . ) ✲ k → k + 1 using f (x) (m k+1 , σ k+1 , C k+1 , . . . ) (αm k , ασ k , C k , . . . ) ✲ k → k + 1 using f (x/α) (αm k+1 , ασ k+1 , C k+1 , . . . ) T (α) ❄ ✻ T -1 (α) T -1 (α) ✻ ❄ T (α)
(m k , σ k , C k , p σ k , p c k ) CMA-ES h (T (m k , σ k , C k , p σ k , p c k )) = T (CMA-ES f (m k , σ k , C k , p σ k , p c k = T -1 (T (m k ,σ k ,C k ,p σ k ,p c k )) )) , (33)
see Fig. 5. Furthermore, for any given f : R n → R, the set of functions H scale ( f )the orbit of f -is an equivalence class.

Proof idea. We investigate the update equations of the state variables comparing the two possible paths from the lower left to the lower right in Fig. 5. The equivalence relation property can be shown elementarily (compare Proposition 1) or using the property that the set {α > 0} is a transformation group over the set {h : R n → R} and therefore induces the equivalence classes H scale ( f ) (see also Proposition 9).

⊓ ⊔

Invariance allows to draw the commutative diagram of Fig. 5. Scale invariance can be interpreted in several ways.

• The choice of scale α is irrelevant for the algorithm, that is, the algorithm has no intrinsic (built-in) notion of scale. • The transformation T in Fig. 5 is a change of coordinate system (here: change of scale) and the update equations are independent of the actually chosen coordinate system, that is, they could be formulated in an algebraic way. • For functions in the equivalence class H scale ( f ), the trace of the algorithm (αm k , ασ k , C k , p σ k , p c k ) will be identical for all k = 0, 1, 2 . . . , given that m 0 and σ 0 are chosen appropriately, for example σ 0 = 1/α and m 0 = σ 0 × a. Then the trace for k = 0 equals (αm 0 , ασ 0 , C 0 , . . . ) = (a, 1, C 0 , . . . ) and the trace does not depend on α for any k ≥ 0.

• From the last point follows that the step-size control has a distinct role in scale invariance. In practice, when α is unknown, adaptation of the step-size that achieves σ k ∝ 1/α can render the algorithm virtually independent of α.

Scale invariance and step-size control also facilitate the possibility of linear convergence in k to the optimum x * , in that

lim k→ ∞ k m k -x * m 0 -x * = exp - c n (34)
exists with c > 0 or equivalently

lim k→ ∞ 1 k log m k -x * = lim k→ ∞ 1 k log m k -x * m 0 -x * = lim k→ ∞ 1 k t ∑ k=1 log m k -x * m k-1 -x * = - c n (35)
and similarly

E log m k+1 -x * m k -x * → - c n for k → ∞ . (36)
Hence, c denotes a convergence rate and for c > 0 the algorithm converges "loglinearly" (in other words, geometrically fast) to the optimum.

In the beginning of this section we had stated two basic demands on a search algorithm, step-size control is meant to address, namely solving linear functions and the sphere function appropriately fast. We now pursue, with a single experiment, whether the demands are satisfied.

Figure 6 shows a run on the objective function f : R n → R, x → x , with n = 20, λ = 12 (the default value, see Table 2) and with σ 0 = 10 -9 chosen far too small given that m 0 = 1. The outcome when repeating this experiment looks always very similar. We discuss the demands in turn.

1. During the first 170 iterations the algorithm virtually "observes" the linear function x → ∑ 20 i=1 x i at point 1 ∈ R 20 . We see during this phase that σ increases geometrically fast (linearly on the log scale). From this observation, and the invariance properties of the algorithm (also rotation invariance, see below), we can safely imply that the demand for linear functions is satisfied. 5), ( 15), ( 25) and ( 28)) on the objective function f : R 20 → R, x → x , as a member of the equivalence class of functions x → g( α xx * ) with identical behavior, given σ 0 ∝ 1/α and m 0 = σ 0 × (const + x * ). Here, m 0 = 1 and the inital step-size σ 0 = 10 -9 is chosen far too small. Left: f (m k ) (think blue graph) and σ k versus iteration number k in a semi-log plot. Right: all components of m k versus k.

m k=600 m k=180 ≈ 10 -9.5 10 0 ≈ exp -c 600 -180 20 . (37)
Solving for c yields c ≈ 1.0 and with min(n, λ ) = λ = 12 we get c ≈ 1.0 ≪ 0.24 = 0.02 min(n, λ ). Our demand on the convergence rate c is more than satisfied. The same can be observed when covariance matrix adaptation is applied additionally (not shown).

The demand on the convergence ( 23) can be rewritten in that

log m k -x * ≈ -c k n + const . ( 38
)
The k in the RHS nominator implies linear convergence in the number of iterations. The n in the denominator implies linear scale-up: the number of iterations to reduce the distance to the optimum by a given factor increases linearly with the dimension n. Linear convergence can also be achieved with covariance matrix adaptation. Given λ ≫ n, linear scale-up cannot be achieved with covariance matrix adaptation alone, because a reliable setting for the learning rate for the covariance matrix is o(1/n). However, step-size control is reliable and achieves linear scale-up given the step-size damping parameter d σ = O(1) in [START_REF] Jägersküpper | Lower bounds for randomized direct search with isotropic sampling[END_REF]. Scale-up experiments are inevitable to support this claim and have been done, for example, in [START_REF] Hansen | Completely derandomized self-adaptation in evolution strategies[END_REF].

Covariance Matrix Adaptation

In the remainder we exploit the f -ranked (i.e. selected and ordered) set (x 1:λ , . . . , x µ:λ ) to update the covariance matrix C. First, we note that the covariance matrix represents variation parameters. Consequently, an apparent principle is to encourage, or reinforce variations that have been successful-just like successful candidate solutions are reinforced in the update of m in [START_REF] Hansen | The CMA evolution strategy: a comparing review[END_REF]. Based on the current set of f -ranked points, the successful variations are (by definition)

x i:λ -m k for i = 1, . . . , µ . (39)
Remark that "successful variation" does not imply f (x i:λ ) < f (m k ) which is neither necessary nor important nor even desirable in general. Even the demand f (x 1:λ ) < f (m k ) would often result in a far too small step-size.

The Rank-µ Matrix

From the successful variations in [START_REF] Teytaud | Lower bounds for evolution strategies using VC-dimension[END_REF] we form a covariance matrix

C µ k+1 = µ ∑ i=1 w i x i:λ -m k σ k × x i:λ -m k T σ k . (40)
Equation ( 40) is analogous to [START_REF] Hansen | The CMA evolution strategy: a comparing review[END_REF] where successful solution points are used to form the new incumbent. We can easily derive the condition

E(C µ k+1 |C k ) = C k (41)
under pure random ranking thus explaining the factors 1/σ k in [START_REF] Teytaud | General lower bounds for evolutionary algorithms[END_REF].

Assuming the weights w i as given, the matrix C µ k+1 maximizes the (weighted) likelihood of the f -ranked steps.

Proposition 7 (Maximum likelihood estimate of C). Given µ ≥ n, the matrix C µ k+1 maximizes the weighted log-likelihood

C µ k+1 = arg max C pos def µ ∑ i=1 w i × log p N x i:λ -m k σ k C , (42)
where p N (x | C) denotes the density of N (0, C) at point x, and therefore the RHS of (42) reads more explicitly

arg max C pos def - 1 2 log det(αC) - 1 2σ k 2 µ ∑ i=1 w i (x i:λ -m k ) T C -1 (x i:λ -m k ) (43)
where α = 2πσ k 2 is irrelevant for the result.

Proof idea. The proof is non-trivial but works similar to the classical non-weighted case.

⊓ ⊔

In contrast to the computation of m in (12), we are not aware of a derivation for optimality of certain weight values in [START_REF] Teytaud | General lower bounds for evolutionary algorithms[END_REF]. Future results might reveal that different weights and/or even a different value for µ are desirable for [START_REF] Hansen | The CMA evolution strategy: A tutorial[END_REF] and [START_REF] Teytaud | General lower bounds for evolutionary algorithms[END_REF]. Before we turn finally to the covariance matrix update, we scrutinize the computation of C µ k+1 .

What is missing?

In Section 3 we have argued to use only the µ best solutions from the last iteration to update distribution parameters. For a covariance matrix update, disregarding the worst solutions might be too conservative and a negative update of the covariance matrix with the µ worst solutions is proposed in [START_REF] Jastrebski | Improving evolution strategies through active covariance matrix adaptation[END_REF]. This idea is not accommodated in this chapter, but has been recently exploited with consistently good results [START_REF] Arnold | Active covariance matrix adaptation for the (1+1)-CMA-ES[END_REF][START_REF] Hansen | Benchmarking a weighted negative covariance matrix update on the bbob-2010 noiseless testbed[END_REF]. An inherent inconsistency with negative updates though is that long steps tend to be worse merely because they are long (and not because they represent a bad direction) meanwhile, unfortunately, long steps also lead to stronger updates.

At first sight we might believe to have covered all variation information given by x i:λm k in the covariance matrix C µ k+1 . On closer inspection we find that the outer product in [START_REF] Teytaud | General lower bounds for evolutionary algorithms[END_REF] removes the sign: using -(x i:λm) instead of x i:λm in (40) yields the same C µ k+1 . One possibility to recover the sign information is to favor the direction x i:λm over -(x i:λm) = m kx i:λ in some way. This seems difficult to accomplish without affecting either the distribution mean (interfering with Proposition 3) or the maximum entropy property. Therefore, we choose a different way to recover the sign information.

Another Evolution Path

We recover the sign information in a classical and rather heuristic way, which turns out to be nevertheless quite effective. We consider an evolution path xm k-j for j > 0, where x might be m k+1 or any x i:λ . We decompose the path into the recent step and the old path

x -m k-j = x -m k + m k -m k-j . ( 44
)
Switching the sign of the last step means using the vector m kx instead of xm k and we get in this case

m k -x + m k -m k-j = 2(m k -x) + x -m k-j = x -m k-j -2(x -m k ) . ( 45
)
Comparing the last line with the LHS of (44), we see that now the sign of the recent step matters. Only in the trivial cases, if either x = m k (zero step) or m k = m k-j (pre-vious zero path) the outer products of ( 44) and ( 45) are identical. Because we will compute the evolution path over a considerable number of iterations j, the specific choice for x should become rather irrelevant and we will use m k+1 in the following.

In practice, we compute the evolution path, analogous to [START_REF] Hansen | Completely derandomized self-adaptation in evolution strategies[END_REF]. We set p c 0 = 0 and use the momentum equation

p c k+1 = (1 -c c ) p c k + h σ c c (2 -c c )µ eff m k+1 -m k σ k , ( 46
)
where

h σ = 1 if p σ k+1 2 < 1 -(1 -c σ ) 2(k+1) 2 + 2
n+1 n and zero otherwise; h σ stalls the update whenever p σ k+1 is large. The implementation of h σ supports the judgment of pursuing a heuristic rather than a first principle here, and is driven by two considerations.

1. Given a fast increase of the step-size (induced by the fact that p σ k+1 is large), the "visible" landscape will change fast and the adaptation of the covariance matrix to the current landscape seems inappropriate, in particular, because 2. the covariance matrix update using p c is asymmetric: a large variance in a single direction can be introduced fast (while p c k+1 is large), but the large variance can only be removed on a significantly longer time scale. For this reason in particular, an unjustified update should be avoided.

While in (46), again, 1c c is the decay factor and 1/c c ≈ (n + 4)/4, the remaining constants are determined by the stationarity condition

p c k+1 ∼ p c k , (47)
given p c k ∼ N (0, C k ) and pure random ranking and h σ ≡ 1. The evolution path p c heavily exploits the sign information. Let us consider, for a given y ∈ R n , two hypothetical situations with m k+1m k = α k y, for k = 0, 1, 2, . . . . We find that for k → ∞

if α k = 1 then p c k → 2 -c c c c y ≈ n + 2 2 y (48) if α k = (-1) k then p c k → (-1) k-1 c c 2 -c c y ≈ (-1) k-1 2 n + 2 y . ( 49
)
Both equations follow from solving the stationarity condition

x = (1 -c c ) × (±x) + c c (2 -c c ) for x.
Combining both equations, we get the ratio between maximal and minimal possible length of p c , given the input vectors have constant length, as

2 -c c c c ≈ n + 2 2 . ( 50
)
Additionally to the matrix C µ k+1 , we use the rank-one matrix p c k+1 p c k+1 T to introduce the missing sign information into the covariance matrix. The update is specified be-low in (51). The update implements the principal heuristic of reinforcing successful variations for variations observed over several iterations.

Evaluation of the Cumulation Heuristic

We evaluate the effect of the evolution path for covariance matrix adaptation. The right plot shows the running lengths for four different parameter settings versus dimension. For n = 3 the smallest speed-up of about 25% is observed for all variants with c -1 c > 1. The speed-up grows to a factor of roughly 2, 4, and 10 for dimensions 10, 30, and 100, respectively, and always exceeds a factor of √ n/2. For c c = 1 (heuristic off) the scaling with the dimension is ≈ n 1.7 . For c -1 c = √ n the scaling becomes ≈ n 1.1 and about linear for c -1 c ≥ n/3. These findings hold for any function, where the predominant task is to acquire the orientation of a constant number of "long axes", in other words to find a few insensitive directions, where yet a large distance needs to be traversed. The assertion in [START_REF] Suttorp | Efficient covariance matrix update for variable metric evolution strategies[END_REF] that c -1 c ∝ n is needed to get a significant scaling improvement turns out to be wrong. For larger population sizes λ , where the rank-µ update becomes more effective, the positive effect reduces and almost vanishes with λ = 10n.

The same experiment has been conducted on other (unimodal) functions. While on many functions the cumulation heuristic is less effective and yields only a rather n-independent and small speed-up (e.g. on the Rosenbrock function somewhat below a factor of two), we have not seen an example yet, where it compromises the performance remarkably. Hence the default choice has become c -1 c ≈ n/4 (see Table 2 in the appendix), because (a) the update for the covariance matrix will have a time constant of c -1 1 ≈ n 2 /2 and we feel that c -1 1 /c -1 c should not be smaller than n, and (b) in our additional experiments the value c -1 c = n is indeed sometimes worse than smaller values.

The Covariance Matrix Update

The final covariance matrix update combines a rank-one update using p c p cT and a rank-µ update using C µ k+1 ,

C k+1 = (1 -c 1 -c µ + c ε ) C k + c 1 p c k+1 p c k+1 T + c µ C µ k+1 , (51)
where p c and C µ k+1 are defined in ( 46) and ( 40) respectively, and

c ε = (1-h σ 2 ) c 1 c c (2- c c
) is of minor relevance and makes up for the loss of variance in case of h σ = 0. The constants c 1 ≈ 2/n 2 and c µ ≈ µ eff /n 2 for µ eff < n 2 are learning rates satisfying c 1 + c µ ≤ 1. The approximate values reflect the rank of the input matrix or the number of input samples, divided by the degrees of freedom of the covariance matrix. The remaining degrees of freedom are covered by the old covariance matrix C k . Again, the equation is governed by a stationarity condition.

Proposition 8 (Stationarity of covariance matrix C). Given pure random ranking and p c k ∼ N (0, C k ) and h σ = 1, we have

E(C k+1 |C k ) = C k . (52)
Proof idea. Compute the expected value of Equation (51).

⊓ ⊔

Finally, we can state general linear invariance for CMA-ES, analogous to scale invariance in Proposition 6 and Fig. 5.

Proposition 9 (Invariance under general linear transformations). The CMA-ES is invariant under full rank linear transformations of the search space, that is, for each f : R n → R invariant under

H GL : f → { f • B -1 : x → f (B -1 x) | B is a full rank n × n matrix} . (53)
The respective bijective state space transformation reads

T B : (m, σ , C, p σ , p c ) → (Bm, σ , BCB T , p σ , Bp c ) . (54)
Furthermore, for each f , the set H GL ( f ) is an equivalence class with identical algorithm trace T B (m k , σ k , C k , p σ k , p c k ) for a state s and the initial state (m 0 , σ 0 , C 0 , p σ 0 , p c 0 ) = T -1 B (s). Proof idea. Straight forward computation of the updated tuple. The equivalence relation property can be shown elementarily (compare Proposition 1) or by recognizing that the set of full rank matrices is a transformation group over the set { f : R n → R} with group action (B, f ) → f • B -1 and therefore induces the equivalence classes H GL ( f ) as orbits of f under the group action.

⊓ ⊔

A commutative diagram, analogous to Fig. 5, applies with T B in place of T (α) and using f (B -1 x) in the lower path. The transformation B can be interpreted as a change of basis and therefore CMA-ES is invariant under linear coordinate system transformations. All further considerations made for scale invariance hold for invariance under general linear transformations likewise.

Because an appropriate (initial) choice of B is usually not available in practice, general linear invariance must be complemented with adaptivity of C to make it useful in practice and eventually adapt a linear encoding [START_REF] Hansen | Adaptive encoding: How to render search coordinate system invariant[END_REF].

Corollary 1 (Adaptive linear encoding and variable metric [START_REF] Hansen | Adaptive encoding: How to render search coordinate system invariant[END_REF]). The covariance matrix adaptation implements an adaptive linear problem encoding, that is, in other words, an adaptive change of basis, or a change of coordinate system, or a variable metric for an evolution strategy.

Proof idea. (The proof can be found in [START_REF] Hansen | Adaptive encoding for optimization[END_REF]). General linear invariance achieves identical performance on f (B -1 x) under respective initial conditions. Here, B is the linear problem encoding used within the algorithm. Changing (or adapting) C without changing m turns out to be equivalent with changing the encoding (or representation) B in a particular way without changing B -1 m (see also [START_REF] Hansen | Invariance, self-adaptation and correlated mutations in evolution strategies[END_REF][START_REF] Hansen | Adaptive encoding for optimization[END_REF]). Also, for each possible encoding we find a respective covariance matrix BB T .

⊓ ⊔

While adaptation of C is essential to implement general linear invariance, rotation invariance does not necessarily depend on an adaptation of C: rotation invariance is already achieved for C ≡ I, because BIB T = I when B is a rotation matrix, compare (54). Nevertheless, it is important to note that covariance matrix adaptation preserves rotation invariance.

Corollary 2 (Rotation invariance). The CMA-ES is invariant under search space rotations.

Proof idea. Rotation invariance follows from Proposition 9 when restricted to the orthogonal group with BB T = I (for any initial state).

⊓ ⊔

An Experiment on Two Noisy Functions

We advocate testing new search algorithms always on pure random, on linear and on various (non-separable) quadratic functions with various initializations. For the (µ/µ W , λ )-CMA-ES this has been done elsewhere with the expected results: parameters are unbiased on pure random functions, the step-size σ grows geometrically fast on linear functions, and on convex quadratic functions the level sets of the search distribution align with the level sets of the cost function, in that C -1 aligns to the Hessian up to a scalar factor and small stochastic fluctuations [START_REF] Hansen | The CMA evolution strategy: a comparing review[END_REF][START_REF] Hansen | Completely derandomized self-adaptation in evolution strategies[END_REF].

Here, we show results on the well-known Rosenbrock function

f (x) = n-1 ∑ i=1 100 (x 2 i -x i+1 ) 2 + (x i -1) 2 ,
where the possible achievement is less obvious. In order to "unsmoothen" the landscape, a noise term is added: each function value is multiplied with

exp α N 2 n × (G +C/10) + α N 2 n × (G +C/10) , (55)
where G and C are standard Gauss (normal) and standard Cauchy distributed random numbers, respectively. All four random numbers in (55) are sampled independently each time f is evaluated. The term is a mixture between the common normal noise 1 + G, which we believe has a principal "design flaw" [31], and the log-normal noise exp(G) which is alone comparatively easy to solve, each mixed with a heavy tail distribution which cannot be alleviated through averaging. We believe that this adds several difficulties on top of each other. We show results for two noise levels, α N = 0.01 and α N = 1. A section through the 5-D and the 20-D landscape for α N = 1 is shown in Fig. 8. The lower dimensional landscape appears more disturbed but is not more difficult to optimize.

Figure 9 shows the output from a typical run for α N = 0.01 of the (µ/µ W , λ )-CMA-ES with m 0 = -1 and σ 0 = 1 (correctly presuming that in all variables m i ± 3σ 0 embrace the optimum at 1). The calling sequence in Matlab was 5opts.evalparallel = 'on'; % only one feval() call per iteration cmaes('frosennoisy', -ones(20,1), 1, opts); % run CMA-ES plotcmaesdat; % plot figures using output files

The default population size for n = 20 is λ = 12. An error of 10 -9 , very close to the global optimum, is reached after about 20 000 function evaluations (without covariance matrix adaptation it takes about 250 000 function evaluations to reach 10 -2 ). The effect of the noise is hardly visible in the performance. In some cases, the T , in some cases the noise leads to a failure to approach any optimum (see also below).

The main challenge on the Rosenbrock function is to follow a winding ridge, in the figure between evaluation 1000 and 15 000. The ridge seems not particularly narrow: the observed axis ratio is about twenty, corresponding to a condition number of 400. But the ridge constantly changes its orientation (witness by the lower right sub-figure). Many stochastic search algorithms are not able to follow this ridge and get stuck with a function value larger than one.

In Fig. 10, the noise term is set to α N = 1 generating a highly rugged landscape (Fig. 8) and making it even harder to follow the winding ridge. Most search algorithms will fail to solve this function 6 . Now, two additional heuristics are examined:

First, restarting the algorithm with increasing population size (IPOP, [START_REF] Auger | A restart CMA evolution strategy with increasing population size[END_REF]). The population size is doubled for each restart. A larger population size λ is more robust to rugged landscapes, mainly because the sample variance can be larger (for µ eff < n, the optimal step-size on the sphere function is proportional to µ eff [START_REF] Arnold | Optimal weighted recombination[END_REF]). Restarting with increasing population sizes is a very effective heuristic when a good termination criterion is available.

Second, applying an uncertainty-handling (UH, [START_REF] Hansen | A method for handling uncertainty in evolutionary optimization with an application to feedback control of combustion[END_REF]). The uncertainty-handling reevaluates a few solutions and measures their resulting rank changes [START_REF] Hansen | A method for handling uncertainty in evolutionary optimization with an application to feedback control of combustion[END_REF]. If the rank changes exceed a threshold, an action is taken. Here, σ is increased. This pre- vents to getting stuck, when the noise disturbs the selection too severely, but it can also lead to divergence. This is of lesser relevance, because in this case the original algorithm would most likely have been stuck anyway. Again, a good termination criterion is essential.

Remark that in both cases, for restarts and with the uncertainty handling, another possible action is to increase the number of function evaluations used for each individual in replacing a single value with a median.

For running IPOP-UH-CMA-ES, the following sequence is added before calling cmaes. Each restart uses the same initial conditions, here m 0 = -1 and σ 0 = 1 from above. For α N = 0.01 (Fig. 9) the uncertainty-handing increases the running length by about 15%, simply due to the reevaluations (not shown). For α N = 1 in Fig. 10, it shortens the running length by a factor of about ten by reducing the number of necessary restarts. Typical for noisy functions, the restart was invoked due to stagnation of the run [START_REF] Hansen | Benchmarking a BI-population CMA-ES on the BBOB-2009 noisy testbed[END_REF]. When repeating this experiment, in about 75% one restart is needed to finally converge to the global optimum with λ = 24. Without uncertaintyhandling it takes usually five to six restarts and a final population size of λ ≥ 384. Without covariance matrix adaptation it takes about 70 times longer to reach a similar precision as in Fig. 10.

Experiments with the well-known Ellipsoid function, In the lower left we can well observe that the algorithm gets stuck "in the middle of nowhere" during the first two launches. See also Fig. 9 for more explanations a population size λ = 48 using the IPOP approach and with the default population size of 12 with UH-CMA-ES.

f (x) = n ∑ i=1 10

Summary

Designing a search algorithm is intricate. We recapitulate the principled design ideas for deriving the CMA-ES algorithm.

• Using a minimal amount of prior assumptions on the cost function f in order to achieve maximal robustness and minimal susceptibility to deceptiveness.

-Generating candidate solutions by sampling a maximum entropy distribution adds the least amount of unwarranted information. This implies the stochastic nature of the algorithm and that no construction of potentially better points is undertaken. This also implies an internal quadratic model-at least when the distribution has finite variances-and stresses the importance of neighborhood. Consequently, a variable neighborhood suggests itself. -Unbiasedness of all algorithm components, given the objective function is random and independent of its argument. This principle suggests that only the current state and the selection information should bias the behavior of the algorithm. Adding another bias would add additional prior assumptions.

We have deliberately violated this principle for uncertainty-handling as used in one experiment, where the step-size is increased under highly perturbed selection. -Only the ranking of the most recently sampled candidate solutions is used as feed-back from the objective function. This implies an attractive invariance property of the algorithm.

Exploiting more specific information on f effectively, for example smoothness, convexity, or (partial) separability will lead to different and more specific design decisions, with a potential advantage on smooth, convex or separable functions respectively.

• Introducing and maintaining invariance properties. Even invariance is related to avoiding prior assumptions as it implies not exploiting specific structure of the objective function f (for example separability). We can differentiate two main cases.

-Unconditional invariance properties do not depend on the initial conditions of the algorithm and strengthen any empirical performance observation. They allow to unconditionally generalize empirical observations to the equivalence class of functions induced by the invariance property. -Invariance properties that depend on state variables of the algorithm (like σ k for scale invariance in Fig. 5) must be complemented with adaptivity. They are particularly attractive, if adaptivity can drive the algorithm fast into the most desirable state. This behavior can be empirically observed for CMA-ES on the equivalence class of convex-quadratic functions.

Step-size control drives step-size σ k close to its optimal value and adaptation of the covariance matrix reduces these functions to the sphere model.

• Exploiting all available information effectively. The available information and its exploitation are highly restricted by the first two demands. Using a deterministic ranking and different weights for updating m and C are due to this design principle. Also the evolution paths in (46) in (51) are governed by exploiting otherwise unused sign information. Using the evolution paths does not violate any of the above demands, but allows to additionally exploit dependency information between successive time steps of the algorithm. • Solving the two most basic continuous domain functions reasonably fast. Solving the linear function and the sphere function reasonably fast implies to introduce step-size control. These two functions are quite opposed: the latter requires convergence, the former requires divergence of the algorithm.

Finally, two heuristic concepts are applied in CMA-ES.

• Reinforcement of the better solutions and the better steps (variations) when updating mean and variances respectively. This seems a rather unavoidable heuristic given a conservative use of information from f . This heuristic bears the maximum likelihood principle. • Orthogonality of successive steps. This heuristic is a rather common conception in continuous domain search.

Pure random search, where the sample distribution remains constant in the iteration sequence, follows most of the above design principles and has some attractive robustness features. However, pure random search neither accumulates information from the past in order to modify the search distribution, nor changes and adapts internal state variables. Adaptivity of state variables however detaches the algorithm from the initial conditions and let (additional) invariance properties come to life. Only invariance to increasing f -value transformations (Proposition 1) is independent of state variables of the search algorithm. We draw the somewhat surprising conclusion that the abstract notion of invariance leads, by advising the introduction of adaptivity, when carefully implemented, to a vastly improved practical performance.

Despite its generic, principled design, the practical performance of CMA-ES turns out to be surprisingly competitive, or even superior, also in comparatively specific problem classes. This holds in particular when more than 100n function evaluations are necessary to find a satisfactory solution [START_REF] Hansen | Comparing results of 31 algorithms from the black-box optimization benchmarking BBOB-2009[END_REF]-even for example on smooth unimodal non-quadratic functions [START_REF] Auger | Experimental comparisons of derivative free optimization algorithms[END_REF], or on highly multimodal functions [START_REF] Hansen | Evaluating the CMA evolution strategy on multimodal test functions[END_REF] and on noisy or highly rugged functions [START_REF] Hansen | Benchmarking a BI-population CMA-ES on the BBOB-2009 noisy testbed[END_REF]. In contrast, much better search heuristics are available given (nearly) convex-quadratic problems or (partially) separable multimodal problems.

Table 1 Summary of the update equations for the state variables in the (µ/µ W , λ )-CMA-ES with iteration index k = 0, 1, 2, . . . . The chosen ordering of equations allows to remove the iteration index in all variables but m k . Unexplained parameters and constants are given in Table 2 Given 2 is symmetric with positive eigenvalues and satisfies

k ∈ N ∪ {0}, m k ∈ R n , σ k ∈ R + , C k ∈ R n×n positive definite, p σ k ∈ R n , and p c k ∈ R n x i ∼ m k + σ k × N (0, C k ) is normally distributed for i = 1, . . . , λ (56)
C k+1 = (1 -c 1 + (1 -h σ 2 ) c 1 c c (2 -c c )) C k + c 1 p c k+1 p c k+1 T + c µ µ ∑ i=1 w i x i:λ -m k σ k × x i:λ -m k T σ k -C k (60)
σ k+1 = σ k × exp 1 ∧ c σ d σ p σ k+1 E N (0, I) -1 (61)
C k -1 2 C k -1

x 2 Fig. 2

 22 Fig. 2 150 samples from a multivariate (standard) normal distribution in 2-D. Both coordinates are i.i.d. according to a standard normal distribution. The circle depicts the one-σ equal density line, the center of the circle is the mean and modal value at zero. In general, lines of equal density (level sets) are ellipsoids. The probability to sample a point outside the dashed box is close to 1 -(1 -2 × 0.0015) 2 ≈ 1/170 x 1

Fig. 5

 5 Fig. 5 Commutative diagram for scale invariance. Vertical arrows depict an invertible transformation (encoding) T of all state variables of CMA-ES with T (α) : (m, σ , C, p σ , p c ) → (αm, ασ , C, p σ , p c ). Horizontal arrows depict one time step of CMA-ES, applied to the respective tuple of state variables. The two possible paths between a state at time k and a state at time k + 1 are equivalent in all (four) cases. For α = 1 the diagram becomes trivial. The diagram suggests that CMA-ES is invariant under the choice of α > 0 in the sense that, given T and T -1 were available, any function x → f (αx) is (at least) as easy to optimization as f

2 . 0 ≈Fig. 6 A

 206 Fig.6A run of CSA-ES (Equations (5), (15), (25) and (28)) on the objective function f : R 20 → R, x → x , as a member of the equivalence class of functions x → g( α xx * ) with identical behavior, given σ 0 ∝ 1/α and m 0 = σ 0 × (const + x * ). Here, m 0 = 1 and the inital step-size σ 0 = 10 -9 is chosen far too small. Left: f (m k ) (think blue graph) and σ k versus iteration number k in a semi-log plot. Right: all components of m k versus k.

Fig. 7

 7 Fig. 7 Number of function evaluations to reach f (x) < 10 -6 on f (x) = x 2 1 + 10 6 ∑ n i=2 x 2 i with m 0 = 1 and σ 0 = 1. For a (backward) time horizon of c -1 c = 1, the cumulation heuristic is, by definition, switched off. Left figure: number of function evaluations, where each point represents a single run, plotted versus the backward time horizon of the evolution path, c -1 c , for n = [3; 10; 30; 100] (from bottom to top). Triangles show averages for c -1 c = √ n and n, also shown on the right. Right figure: average number of function evaluations divided by n, from [10; 3; 2; 1] = ⌊10/⌊ √ n⌋⌋ runs, plotted versus n for (from top to bottom) c -1 c = 1; √ n; n+3 3 ; n. Compared to c c = 1, the speed-up exceeds in all cases a factor of √ n/2 (dashed line).

Fig. 8

 8 Fig. 8 Both figures show three sections of the Rosenbrock function for α N = 1 and argument x = β × 1 + 1 20 N (0, I). All graphs show 201 points for β ∈ [-0.5, 1.5] and a single realization of N (0, I) in each subfigure. The left basin about zero is initially highly attractive (compare e.g. Fig. 9, upper right) but is not nearby a local or global optimum. The basin around β = 1 is close to the global optimum at 1 and monotonically (non visibly) connected to the left basin

  (f), cyan:f-min(f), green:sigma, red:

Fig. 9 A

 9 Fig. 9 A typical run of the (µ/µ W , λ )-CMA-ES on the Rosenbrock function (n = 20) with a small disturbance of the function value (α N = 0.01). All values are plotted against number of objective function evaluations. Upper left: iteration-wise best function value (thick blue graph), median and worst function value (black graphs, mainly hidden), square root of the condition number of C k (increasing red graph), smallest and largest coordinate-wise standard deviation of the distribution N 0, σ k 2 C k with final values annotated (magenta), and σ k following closely the largest standard deviation (light green). Lower left: square roots of eigenvalues of C k , sorted. Upper right: incumbent solution m k . Lower right: square roots of diagonal elements of C k

Fig. 10 A

 10 Fig. 10 A typical run of the IPOP-UH-CMA-ES on the noisy Rosenbrock function (n = 20, α N = 1), a (µ/µ W )-CMA-ES with Uncertainly Handling restarted with Increasing POPulation size. The highly rugged lines, partly beyond 10 5 , in the upper left depict the worst measured function value (out of λ ). One restart was necessary to converge close to the global optimum. See also Fig. 9 for more explanations

6 i- 1 n-1 x 2 iFig. 11

 1211 Fig. 11 Two typical runs of the IPOP-CMA-ES (left) and UH-CMA-ES (right, with uncertaintyhandling) on the noisy Ellipsoid function (n = 20, α N = 1). With α N = 0 the Ellipsoid is solved in about 22000 function evaluations.In the lower left we can well observe that the algorithm gets stuck "in the middle of nowhere" during the first two launches. See also Fig.9for more explanations

m 1 2m

 1 k+1 = m k + c m µ ∑ i=1 w i (x i:λm k ) where f (x 1:λ ) ≤ • • • ≤ f (x µ:λ ) ≤ f (x µ+1:λ ) . . . (57)p σ k+1 = (1c σ ) p σ k + c σ (2c σ )µ eff C k -k+1m k c m σ k(58)p c k+1 = (1c c ) p c k + h σ c c (2c c )µ eff m k+1m k c m σ k (59)

  c σ ) 2(k+1) < 2 + 4 n+1 n 0 otherwise , for E N (0, I) = √ 2 Γ( n+1 2 )/Γ( n 2 ) ≈ n -1/2 we use the better approximation √ n 1 -1 4n + 1 21n 2 , and C k -1

3

  Commutative diagram for invariance. Vertical arrows depict an invertible transformation (encoding) T of the state variables. Horizontal arrows depict one time step of algorithm A , using the respective function and state variables. The two possible paths between a state at time k and a state at time k + 1 are equivalent in all (four) cases. The two paths from upper left to lower right are reflected in Equation 8. For f = h the diagram becomes trivial with T f →h as the identity. One interpretation of the diagram is that given T -1f →h , any function h can be optimized like f

That is, to find a sequence θ k , k = 1,

2, 3 . . . , such that lim k→∞ E( f (x|θ k )) = f * .

Different learning rates might be related to some parameters in the distribution being orthogonal.

Source code is available at http://www.lri.fr/˜hansen/cmaes inmatlab.html and will be accessible at http://cma.gforge.inria.fr/ in future. In our experiment, version 3.40.beta was used with Matlab.

There is a simple way to smoothen the landscape: a single evaluation can be replaced by the median (not the mean) of a number of evaluations. Only a few evaluations reduce the dispersion considerably, but about 1000 evaluations are necessary to render the landscape similarly smooth as with α N = 0.01. Together with (µ/µ W )-CMA-ES, single evaluations, as in Fig.10, need overall the least number of function evaluations (comprising restarts).

= (C k ) -1 . The binary ∧ operator depicts the minimum of two values with low operator precedence. The default parameter values are shown in Table2.

Acknowledgements The authors would like to express their gratitude to Marc Schoenauer for his kind and consistent support.

Appendix

The (µ/µ W , λ )-CMA-ES, as described in this chapter, is summarized in Table 1. We have p σ k=0 = p c k=0 = 0, C k=0 = I, while m k=0 ∈ R n and σ k=0 > 0 are user defined. Additionally, x i:λ is the i-th best of the solutions x 1 , . . . , x λ , Nikolaus Hansen and Anne Auger

⌊3 ln n⌋ population size, see also [START_REF] Auger | A restart CMA evolution strategy with increasing population size[END_REF][START_REF] Hansen | Benchmarking a BI-population CMA-ES on the BBOB-2009 function testbed[END_REF] µ = λ 2 parent number

cumulation constant for step-size, 1/c σ is the respective time constant

covariance matrix learning rate for the rank one update using p c

covariance matrix learning rate for rank-µ update α cov = 2 ∧ λ /3 could be chosen < 2, e.g. α cov = 0.5 for noisy problems


Emmanuel Castanier

Clement Jonquet

Soumia Melzi

Pierre Larmande

Manuel Ruiz

Patrick Valduriez

Semantic Annotation Workflow using Bio-Ontologies

niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Semantic Annotation Workflow using Bio-Ontologies

• Second, annotations are semantically expanded using knowledge of the ontologies.

• Third, all annotations are aggregated and scored according to the context in which they have been created.

Introduction

Biologists


Julien Grosjean

Lina F Soualmia

Khedidja Bouarech

Clément Jonquet

Stéfan J Darmoni

An Approach to Compare Bio-Ontologies Portals

Keywords: Controlled Vocabulary, Internet, Information Storage and Retrieval, Terminology as subject, Health terminologies and ontologies, Crosslingual access to health information, Portal

Background: main biomedical information retrieval systems are based on controlled vocabularies and most specifically on terminologies or ontologies (T/O). These classification structures allow indexing, coding, annotating different kind of documents. Many T/O have been created for different purposes and it became a problem for finding specific concepts in the multitude of existing nomenclatures. The NCBO (National Center for Biomedical Ontologies) BioPortal 2 and the CISMeF (Catalogue et Index des Sites Médicaux de langue Franc ¸aise) HeTOP 3 projects have been developed to tackle this issue. Objective: the present work consists in comparing both portals. Methods: we hereby are proposing a set of criteria to compare bio-ontologies portals in terms of goals, features, technologies and usability. Results: BioPortal and HeTOP have been compared based on the given criteira. While both portals are designed to store and make T/O available to the community and are sharing many basic features, they are also very different mainly because of their basic purposes. Conclusion: thanks to the comparison criteria, we can assume that a merge between BioPortal and HeTOP is possible in terms of functionnalities. The main difficulties will be about merging the data repositories and applying different policies on T/O content.

Introduction

Because of the constant growth of biomedical data it became mandatory to index or annotate it with controlled and structured vocabularies in order to store it and to retrieve it intelligently. A key aspect in addressing semantic interoperability for biomedical data is the use of terminologies or ontologies (T/O) as a common data structure data [START_REF] Bodenreider | Bio-ontologies: Current Trends and Future Directions[END_REF] [START_REF] Rubin | Biomedical ontologies: a functional perspective[END_REF]. Many T/O have been created in the last decade for different purposes: indexing or annotating documents, organising knowledge, inferring facts, etc. Several tools have been created to store, search and use multiple T/O simultaneously. Among them, the UMLS (Unified Medical Language System) [START_REF] Bodenreider | The Unified Medical Language System (UMLS): integrating biomedical terminology[END_REF], the EBI Ontology Lookup Service [START_REF] Côté | The Ontology Lookup Service: bigger and better[END_REF], the NCBO BioPortal [START_REF] Rubin | National Center for Biomedical Ontology: Advancing Biomedicine through Structured Organization of Scientific Knowledge[END_REF] or the CISMeF HeTOP [START_REF] Grosjean | Health Multi-Terminology Portal: a semantics added-value for patient safety. Patient Safety Informatics -Adverse Drug Events, Human Factors and IT Tools for Patient Medication Safety[END_REF].

The NCBO project (National Center for Biomedical Ontology) of the university of Stanford and the CISMeF (Catalogue et Index des Sites Médicaux de langue Franc ¸aise, Rouen University Hospital, France) have invested considerable efforts in the development of tools and services based on T/O to assist health professionals to search resources on the Internet and to use T/O. Both groups have developed web portals, respectively BioPortal and HeTOP that offer various services to find, index, browse, visualise and annotate T/O. This article aims to compare BioPortal and HeTOP by defining objective comparison points based on their functionalities and features. Thus, it could be possible to define a strategy to merge both portals into one unique solution, offering the best services for human users and programmes.

Material & Methods

BioPortal -http://bioportal.bioontology.org

Developed by the NCBO, BioPortal is a repository of biomedical ontologies which hosts more than 350 ontologies in different formats [START_REF] Noy | BioPortal: ontologies and integrated data resources at the click of a mouse[END_REF], [START_REF] Whetzel | Bio-Portal: enhanced functionality via new Web services from the National Center for Biomedical Ontology to access and use ontologies in software applications[END_REF]. These ontologies are regularly updated by users and accessible via a web site for the humans and web-based services for programmes. BioPortal is a library of community ontologies [START_REF] Noy | Where to Publish and Find Ontologies? A Survey of Ontology Libraries[END_REF] designed as a "onestop shop" repository. Users have access to the ontologies with or without restriction (depending on the level of restriction set by the editor) and can access to editing, comment, rating and content adding operations.

HeTOP -http://www.hetop.org

HeTOP (Health Terminology/Ontology Portal) [START_REF] Grosjean | Health Multi-Terminology Portal: a semantics added-value for patient safety. Patient Safety Informatics -Adverse Drug Events, Human Factors and IT Tools for Patient Medication Safety[END_REF], is a T/O portal developed by the CIS-MeF team. It hosts more than 55 T/O in several languages. Most of the T/O are international or French national references such as MeSH, ICD-10, or CCAM. These T/O are regularly updated and are accessible via a web site and a web-based service. HeTOP has been designed as reference multi-terminology and multi-lingual portal [START_REF] Grosjean | Multi-terminology cross-lingual model to create the Health Terminology/Ontology Portal[END_REF] to help librarians, translators, students and health professionals to retrieve resources and knowledge across a high variety of complex health fields.

Criteria and Comparison Approach

The work of [START_REF] Noy | Where to Publish and Find Ontologies? A Survey of Ontology Libraries[END_REF] has inspired us for establishing a criteria list for comparing BioPortal and HeTOP properties (technologies, methodologies, policies, target users, final purposes, features,. . . ) but also from other similar portals. We categorised all the criteria in four groups: (i) Content, (ii) Functions & Tools, (iii) User Interface & Usability, (iv) Methods & Technologies.

Results

Content Comparison

The comparison of content reveals 3 important differences between BioPortal and HeTOP: (i) The volume of data is more important in BioPortal considering the T/O num-bers or concept numbers (respectively 5,960,457 and 1,951,834). However, since T/O have various number of concepts and terms (terms are preferred terms plus synonyms, in any language) and since HeTOP is dealing with multilingual content, T/O numbers are not suited comparison indicators: the more relevant figures are the number of terms and the total number of relations. Indeed, BioPortal and HeTOP have around the same number of terms (about 6,600,000). Unfortunately, it is not possible to easily calculate the total number of relations in BioPortal. (ii) T/O formats and update frequencies are quite different in BioPortal compared to HeTOP. While BioPortal is focused on ontologies and takes advantage of standard formats and programmes, HeTOP is also hosting heterogeneous representation formats such as Microsoft Excel files, XML files or database dumps. A special work has to be done for every new T/O source: this can not be executed automatically and it implies expertise and development. (iii) One major difference between BioPortal and HeTOP is the expertise brought to every T/O. While BioPortal is automatically importing ontologies and does not change anything, except new automatic mappings and some manual users mappings, each HeTOP hosted T/O undergoes a series of process to leverage its content and meta-data (new translations, synonyms, definitions, relations, mappings, etc.).

Functions & Tools Comparison

Despite the similarity of basic tools for both portals, some details and tools differ slightly but have a direct consequence for the human users. (i) The BioPortal search engine only searches for exact terms in English (among preferred terms and synonyms) while HeTOP search engine is able to add wildcards to search for terms containing the query, in two languages simultaneously. This has a direct impact on how users can search terms; for instance, if one searches "myopathy" (in English) in the NCIT 4 , BioPortal retrieves 5 terms whereas HeTOP 25; because wildcards are managed, HeTOP search engine is actually querying "*myopathy*" and retrieves terms such as "Cardiomyopathy".

User Interface & Usability Comparison

To compare the user experience using both portals, we compared two functions involving time responses. First, a comparison have been made between search engines performances in terms of response time and result numbers; we picked and performed 10 random queries on both search engines and we measured the user experience time (with the FireBug Mozilla Firefox plug-in) and noted the results. No options have been selected in both portals and the searches have been made among all T/O, in English and in French in HeTOP and only in English in BioPortal (no possible multilingual search). Results are 5.57 sec/19.7 (average response time/number of results) and 3.17 sec/359.4 respectively for BioPortal and HeTOP. About multilingualism, within BioPortal, T/O in other languages than English are mostly available as "views" of the corresponding English T/O (e.g., the French MeSH is a view of the MeSH) but it is impossible to get the French term while browsing the English one and vice-versa. Within HeTOP, T/O are not language specific (e.g., the MeSH exists only once with the available translated terms). Therefore it is easy from the English term to get to the French one and vice-versa. The whole user interface is internationalised and the search can be performed per language. Switching from one language to another is context sensitive.

Methods & Technical Comparison

While BioPortal is based on a RDF triple store data model, HeTOP has a metamodel for T/O which encapsulates specific T/O models into an Oracle 11g r2 relational database. Both portals are coupled to a web-based service (http://www.bioontology. org/wiki/index.php/Resource_Index_REST_Web_Service_User_Guide, http: //cispro.chu-rouen.fr/CISMeFhetopservice/. About code license, APIs are open for BioPortal and proprietary for HeTOP.

Conclusion

As described above, BioPortal and HeTOP are sharing many features such as T/O browsing, tab representation and resources access tools. Both portals are valuable tools to support research projects. However, their policies and basic purposes are significantly different. While BioPortal is opened to the community T/O, HeTOP is focusing on reference T/O with experts interventions. Thus, BioPortal only accepts ontologies (at least files in ontology formats) whereas HeTOP provides access to non-ontological sources. Bio-Portal has been created as a "one-stop shop" repository and can be instantiated in other environments using the virtual appliance. On the contrary, HeTOP is a static server/web site acting such as a platform to help various kinds of users for different goals. Moreover, HeTOP is focusing on T/O content, working hard with experts and curators to leverage lexicons, mappings and other knowledge resources. This has a direct consequence on data volume and restriction policies related to download, versioning or updates. Several criteria of this comparison study are focusing on usability: a portal dedicated to be used by experts or lay people has to be understandable and usable with less efforts and knowledge. A special work is made on HeTOP to deal with it: many T/O meta-data labels (attributes, relations, ...) are translated in several languages and well defined. On the other hand, BioPortal is more open than HeTOP: users can upload ontologies and annotate concepts and they have more options. BioPortal proposes persistent concept URL thanks to the REST technology. However, HeTOP is more user-oriented with high quality and multilingual T/O content, it is also faster (search engine and navigation) and maybe more adapted to lay people and students which means that it is more reliable and more friendly to use on a daily basis. Unfortunately, no evaluation on user's satisfaction has been performed on HeTOP nor BioPortal; it is an ongoing work for HeTOP with a set of online questions. This study helps to understand both portals philosophies and functionalities. Furthermore, it allows to point at BioPortal and HeTOP advantages and drawbacks. Despite the differences, we can assume that a merge is possible. The vast majority of policy points and functionalities could be kept in both approaches to create a single portal. The main difficulty would be the technical choices to merge a RDF data store and a relational database with many specificities. To tackle this, it would be possible to integrate the two data layers in a unique system coupled to a single API. This merge would have a great cost but it would be a considerable benefit for the biomedical community and T/O users. This work is part of the SIFR (Semantic Indexing of French Biomedical Data Resources) research project 5 , in collaboration betweem LIRMM, CISMeF and NCBO.

Corresponding author: julien.grosjean@chu-rouen.fr

National Cancer Institute thesaurus, edited by the National Cancer Institute: http://ncit.nci.nih. gov/

http://www.lirmm.fr/sifr/

Acknowledgements

This work is issued from a collaboration between NCBO, CISMeF and the LIRMM. We used the knowledge of the partners, the scientific productions and the technical documentations of BioPortal and HeTOP to define comparison criteria. We thank the NCBO team for its cooperation on technical elements. This work was supported in part by the French National Research Agency under JCJC program, grant ANR-12-JS02-01001, as well as by University Montpellier 2, CNRS and IBC project.