The (1+$\lambda$) Evolutionary Algorithm with Self-Adjusting Mutation Rate

TL;DR

This paper introduces a simple self-adjusting mutation rate mechanism for the $(1+\lambda)$ evolutionary algorithm, demonstrating improved optimization efficiency on the OneMax problem and potential for broader application.

Contribution

It presents a novel, parameter-free self-adjusting mutation rate method that outperforms traditional fixed-rate approaches in evolutionary algorithms.

Findings

Expected optimization time is $O(n\lambda/\log\lambda + n\log n)$

The method outperforms classic $(1+\lambda)$ EA in efficiency

Performance matches the best possible among all $\lambda$-parallel black-box algorithms

Abstract

We propose a new way to self-adjust the mutation rate in population-based evolutionary algorithms in discrete search spaces. Roughly speaking, it consists of creating half the offspring with a mutation rate that is twice the current mutation rate and the other half with half the current rate. The mutation rate is then updated to the rate used in that subpopulation which contains the best offspring. We analyze how the $(1+\lambda)$ evolutionary algorithm with this self-adjusting mutation rate optimizes the OneMax test function. We prove that this dynamic version of the $(1+\lambda)$ EA finds the optimum in an expected optimization time (number of fitness evaluations) of $O(n\lambda/\log\lambda+n\log n)$. This time is asymptotically smaller than the optimization time of the classic $(1+\lambda)$ EA. Previous work shows that this performance is best-possible among all $\lambda$-parallel mutation-based unbiased black-box algorithms. This result shows that the new way of adjusting the mutation rate can find optimal dynamic parameter values on the fly. Since our adjustment mechanism is simpler than the ones previously used for adjusting the mutation rate and does not have parameters itself, we are optimistic that it will find other applications.

Equations

$$\text{textsc}OneMax:\{0,1{\}}^n\to \mathbb{R};(x_1,\dots ,x_n)\mapsto \sum _{i=1}^n{x}_i.$$

$$E(T)=\Theta \left(\frac{n}{\log \lambda }+\frac{n\log n}{\lambda }\right).$$

$$E(T)=O\left(\frac{n}{\log \lambda }+\frac{n\log n}{\sqrt{\lambda }}\right).$$

$$E(X_t-X_{t+1}\mid X_t=x)\ge h(x)$$

$$E(T\mid X_0=x^{\prime })\le \frac{x_{\min }}{h(x_{\min })}+{\int }_{x_{\min }}^{x^{\prime }}\frac{1}{h(x)}\mathrm{d}x.$$

$$E(X_t-X_{t+1}\mid X_t=x)\ge \delta x,$$

$$E(T\mid X_0=x^{\prime })\le \frac{1+\ln \left(\frac{x^{\prime }}{x_{\min }}\right)}{\delta }.$$

$$\mathrm{Pr}(X\ge E(X)+\Delta )$$

$$\mathrm{Pr}(X_t\ge b)\le 2e^{\frac{2d}{3c}(1-b/c)}.$$

$$\Omega\bigg{(}\frac{\log\lambda}{\log\frac{en}{k}}\bigg{)}$$

$$\bigg{[}\frac{1}{2\ln\frac{en}{k}}\,\ln\lambda,\frac{4n^{2}}{(n-2k)^{2}}\,\ln\lambda\bigg{]}$$

$$q(k,i,r)=\sum _{j=0}^{k-i}\left(\genfrac{}{}{0ex}{}{k}{i+j}\right)\left(\genfrac{}{}{0ex}{}{n-k}{j}\right){\left(\frac{r}{n}\right)}^{i+2j}{\left(1-\frac{r}{n}\right)}^{n-i-2j}$$

$$Q(k,2r,2r)\ge q(k,2r,2r)\ge \left(\genfrac{}{}{0ex}{}{k}{2r}\right){\left(\frac{2r}{n}\right)}^{2r}{\left(1-\frac{2r}{n}\right)}^n$$

$$\ge 4{\left(\frac{k}{2r}\right)}^{2r}{\left(\frac{2r}{n}\right)}^{2r}{e}^{-2r}\ge 4{\left(\frac{k}{en}\right)}^{2r}\ge 4{\left(\frac{k}{en}\right)}^{2c_1(k)\ln \lambda }=\frac{4}{\lambda }.$$

$$\frac{q(k,i,2r)}{q(k,i,r/2)}\geq 4^{i}\left(\frac{1-2r/n}{1-r/(2n)}\right)^{n}\geq\big{(}1-o(1)\big{)}4^{i}e^{-1.5r}.$$

$$\frac{q(k,i,2r)}{q(k,i,r/2)}\geq\big{(}1-o(1)\big{)}\left(\frac{4}{e}\right)^{2r}>4.$$

$$(1-0.2)\cdot 4/6=\mathrm{0.64.}$$

$$E(X(k,r))=kp-(n-k)p=\frac{(2k-n)r}{n},$$

$$\mathrm{Var}(X(k,r))=np(1-p)=r\left(1-\frac{r}{n}\right).$$

$$\mathrm{Pr}(X\ge E(X)+\Delta )\le \exp \left(\frac{-{\Delta }^2}{2\mathrm{Var}(X)+2\Delta /3}\right)\le \exp \left(\frac{-{\Delta }^2}{2\mathrm{Var}(X)+2\Delta }\right)$$

$$\mathrm{Pr}(X(k,2r)\ge E(X(k,r/2)))\le \exp \left(\frac{-{\Delta }^2}{2\mathrm{Var}(X(k,2r))+2\Delta }\right)$$

$$=\exp \left(\frac{-9(n-2k{)}^{2}r}{4n(7n-8r-6k)}\right)\le \exp \left(-\frac{9(n-2k{)}^2{c}_2(k)\ln \lambda }{28n^2}\right)=\frac{1}{{\lambda }^{9/7}}.$$

$$X^{+}$$

$$X^{-}$$

$$\mathrm{Pr}(X(k,r/2)\ge E(X(k,r/2)))$$

$$\ge 0.037\cdot 0.5=\mathrm{0.0185.}$$

$$\mathrm{Pr}(X(k,r)\ge 0)$$

$$=\exp \left(-\frac{(n-2k{)}^{2}r}{2n(2n-2k-r)}\right)\le \exp \left(-\frac{(n-2k{)}^{2}r}{4n^2}\right).$$

$$B(x):=\mathrm{Pr}(X^{+}=x)=\left(\genfrac{}{}{0ex}{}{k}{x}\right)p^x(1-p{)}^{k-x}\text{ for all }x\in \{0,1,\dots ,k\},$$

$$F(y):=\mathrm{Pr}(X^{+}\ge y)=\sum _{i=\lceil y\rceil }^{k}B(i)\text{ for all }y\in [0,k].$$

$$\frac{B(x+1)}{B(x)}=\frac{k-x}{x+1}\cdot \frac{p}{1-p}\le \frac{u-up}{u-up+1-p}<1,$$

$$\frac{B(\lceil u+2\delta \rceil )}{B(\lceil u\rceil )}$$

$$\ge {\left(\frac{k-u-2\delta }{k(1-p)}\right)}^{2\delta }\frac{u^{2\delta }}{\lceil u+1\rceil \cdots \lceil u+2\delta \rceil }\ge {\left(1-\frac{2\delta }{3u}\right)}^{2\delta }\frac{u^{2\delta }}{\lceil u+1\rceil \cdots \lceil u+2\delta \rceil },$$

$$\lceil u+2\delta \rceil !\le \sqrt{2\pi \lceil u+2\delta \rceil }{\left(\frac{\lceil u+2\delta \rceil }{e}\right)}^{\lceil u+2\delta \rceil }\exp \left(\frac{1}{12(\lceil u+2\delta \rceil )}\right),$$

Full text

We analyze how the $(1+\lambda)$ evolutionary algorithm with this self-adjusting mutation rate optimizes the OneMax test function. We prove that this dynamic version of the $(1+\lambda)$ EA finds the optimum in an expected optimization time (number of fitness evaluations) of $O(n\lambda/\log\lambda+n\log n)$. This time is asymptotically smaller than the optimization time of the classic $(1+\lambda)$ EA. Previous work shows that this performance is best-possible among all $\lambda$-parallel mutation-based unbiased black-box algorithms.

This result shows that the new way of adjusting the mutation rate can find optimal dynamic parameter values on the fly. Since our adjustment mechanism is simpler than the ones previously used for adjusting the mutation rate and does not have parameters itself, we are optimistic that it will find other applications.

1 Introduction

Evolutionary algorithms (EAs) have shown a remarkable performance in a broad range of applications. However, it has often been observed that this performance depends crucially on the use of the right parameter settings. Parameter optimization and parameter control are therefore key topics in EA research. Since these have very different characteristics in discrete and continuous search spaces, we discuss in this work only evolutionary algorithms for discrete search spaces.

Theoretical research has contributed to our understanding of these algorithms with mathematically founded runtime analyses, many of which show how the runtime of an EA is determined by its parameters. The majority of these works investigate static parameter settings, i. e., the parameters are fixed before the start of the algorithm and are not changed during its execution. More recently, a number of results were shown which prove an advantage of dynamic parameter settings, that is, the parameters of the algorithm are changed during its execution. Many of these rely on making the parameters functionally dependent on the current state of the search process, e.g., on the fitness of the current-best individual. While this provably can lead to better performances, it leaves the algorithm designer with an even greater parameter setting task, namely inventing a suitable functional dependence instead of fixing numerical values for the parameters. This problem has been solved by theoretical means for a small number of easy benchmark problems, but it is highly unclear how to find such functional relations in the general case.

A more designer-friendly way to work with dynamic parameters is to modify the parameters based on simple rules taking into account the recent performance. A number of recent results shows that such on the fly or self-adjusting parameter settings can give an equally good performance as the optimal fitness-dependent parameter setting, however, with much less input from the algorithm designer. For example, good results have been obtained by increasing or decreasing a parameter depending on whether the current iteration improved the best-so-far solution or not, e.g., in a way resembling the $1/5$-th rule from continuous optimization.

Such success-based self-adjusting parameter settings can work well when there is a simple monotonic relation between success and parameter value, e.g., when one speculates that increasing the size of the population in an EA helps when no progress was made. For parameters like the mutation rate, it is not clear what a success-based rule can look like, since a low success rate can either stem from a too small mutation rate (regenerating the parent with high probability) or a destructive too high mutation rate. In [19], a relatively complicated learning mechanism was presented that tries to learn the right mutation strength by computing a time-discounted average of the past performance stemming from different parameter values. This learning mechanism needed a careful trade-off between exploiting the currently most profitably mutation strength and experimenting with other parameter values and a careful choice of the parameter controlling by how much older experience is taken less into account than more recent observations.

1.1 A New Self-Adjusting Mechanism for Population-Based EAs

In this work, we propose an alternative way to adjust the mutation rate on the fly for algorithms using larger offspring populations. It aims at overcoming some of the difficulties of the learning mechanism just described. The simple idea is to create half the offspring with twice the current mutation rate and the other half using half the current rate. The mutation rate is then modified to the rate which was used to create the best of these offspring (choosing the winning offspring randomly among all best in case of ambiguity). We do not allow the mutation rate to leave the interval $[2/n,1/4]$, so that the rates used in the subpopulations are always in the interval $[1/n,1/2]$.

We add one modification to the very basic idea described in the first paragraph of this section. Instead of always modifying the mutation rate to the rate of the best offspring, we shall take this winner’s rate only with probability a half and else modify the mutation rate to a random one of the two possible values (twice and half the current rate). Our motivation for this modification is that we feel that the additional random noise will not prevent the algorithm from adjusting the mutation rate into a direction that is more profitable. However, the increased amount of randomness may allow the algorithm to leave a possible basin of attraction of a locally optimal mutation rate. Observe that with probability $\Theta(1/n^{2})$, a sequence of $\log_{2}n$ random modifications all in the same direction appears. Hence with this inverse-polynomial rate, the algorithm can jump from any mutation rate to any other (with the restriction that only a discrete set of mutation rates can appear). We note that the existence of random modifications is also exploited in our runtime analysis, which will show that the new self-adjusting mechanism selects mutation rates good enough to lead to the asymptotically optimal runtime among all dynamic choices of the mutation rate for the (1+$\lambda$) EA.

In this first work proposing this mechanism, we shall not spend much effort fine-tuning it, but rather show in a proof-of-concept manner that it can find very good mutation rates. In a real application, it is likely that better results are obtained by working with three subpopulations, namely an additional one using (that is, exploiting) the current mutation rate. Also, it seems natural that more modest adjustments of the mutation rate, that is, multiplying and dividing the rate by a number $F$ that is smaller than the value $F=2$ used by our mechanism, is profitable. We conduct some elementary experiments supporting this intuition in Section 9.

1.2 Runtime Analysis for the Self-Adjusting (1+$\lambda$) EA on OneMax

To prove that the self-adjusting mechanism just presented can indeed find good dynamic mutation rates, we conduct a rigorous runtime analysis for our self-adjusting (1+$\lambda$) EA on the classic test function

[TABLE]

The runtime of the (1+$\lambda$) EA with fixed mutation rates on OneMax is well understood [13, 27]. In particular, Gießen and Witt [27] show that the expected runtime (number of generations) is $(1\pm o(1))\left(\frac{1}{2}\cdot\frac{n\ln\ln\lambda}{\ln\lambda}+\frac{e^{r}}{r}\cdot\frac{n\ln n}{\lambda}\right)$ when a mutation rate of $r/n$, $r$ a constant, is used. Thus for $\lambda$ not too large, the mutation rate determines the leading constant of the runtime, and a rate of $1/n$ gives the asymptotically best runtime.

As a consequence of their work on parallel black-box complexities, Badkobeh, Lehre, and Sudholt [2] showed that the (1+$\lambda$) EA with a suitable fitness-dependent mutation rate finds the optimum of OneMax in an asymptotically better runtime of $O(\frac{n}{\log\lambda}+\frac{n\log n}{\lambda})$, where the improvement is by a factor of $\Theta(\log\log\lambda)$.111As usual in the analysis of algorithms, we write $\log$ when there is no need to specify the base of the logarithm, e.g., in asymptotic terms like $\Theta(\log n)$. To avoid the possible risk of an ambiguity for small arguments, as usual, we set $\log x=\max\{1,\log x\}$. All other logarithms are used in their classic meaning, that is, for all $x>0$, we denote by $\ln x$ the natural logarithm of $x$ and by $\log_{2}x$ its binary logarithm. This runtime is best-possible among all $\lambda$-parallel unary unbiased black-box optimization algorithms. In particular, no other dynamic choice of the mutation rate in the (1+$\lambda$) EA can achieve an asymptotically better runtime. The way how the mutation rate depends on the fitness in the above result, however, is not trivial. When the parent individual has fitness distance $d$, then mutation rate employed is $p=\max\{\frac{\ln\lambda}{n\ln(en/d)},\frac{1}{n}\}$.

Our main technical result is that the (1+$\lambda$) EA adjusting the mutation rate according to the mechanism described above has the same (optimal) asymptotic runtime. Consequently, the self-adjusting mechanism is able to find on the fly a mutation rate that is sufficiently close to the one proposed in [2] to achieve asymptotically the same expected runtime.

Theorem 1**.**

Let $\lambda\geq 45$ and $\lambda=n^{O(1)}$. Let $T$ denote the number of generations of the (1+$\lambda$) EA with self-adjusting mutation rate on OneMax. Then,

[TABLE]

This corresponds to an expected number of function evaluations of $\Theta(\frac{\lambda n}{\log\lambda}+n\log n)$.

To the best of our knowledge, this is the first time that a simple mutation-based EA achieves a super-constant speed-up via a self-adjusting choice of the mutation rate.

As an interesting side remark, our proofs reveal that a quite non-standard but fixed mutation rate of $r=\ln(\lambda)/2$ also achieves the $\Theta(\log\log\lambda)$ improvement as it implies the bound of $\Theta(n/\log\lambda)$ generations if $\lambda$ is not too small. Hence, the constant choice $r=O(1)$ as studied in [27] does not yield the asymptotically optimal number of generations unless $\lambda$ is so small that the $n\log n$-term dominates.

Lemma 2**.**

Let $\lambda\geq 45$ and $\lambda=n^{O(1)}$, Let $T$ denote the number of generations of the (1+$\lambda$) EA with fixed mutation rate $r=\ln(\lambda)/2$. Then,

[TABLE]

This corresponds to an expected number of function evaluations of $O(\frac{\lambda n}{\log\lambda}+\sqrt{\lambda}n\log n)$.

The paper is structured as follows: In Section 2 we give an overview over previous analyses of the (1+$\lambda$) EA and of self-adjusting parameter control mechanism in EAs from a theoretical perspective. In Section 3 we give the algorithm and the mutation scheme. For convenience, we also state some key theorems that we will frequently use in the rest of the paper. Afterwards, Section 4 gives a high-level overview of our main proof strategy. The following technical sections deal with the runtime analysis of the expected time spent by the (1+$\lambda$) EA on OneMax in each of three regions of the fitness distance $d$. We label these regions the far region, middle region and near region, each of which will be dealt with in a separate section. The proof of the main theorem and of Lemma 2 is then given in Section 8. Finally, we conclude in Section 10.

2 Related Work

Since this is a theoretically oriented work on how a dynamic parameter choice speeds up the runtime of the (1+$\lambda$) EA on the test function OneMax, let us briefly review what is known about the theory of this EA and dynamic parameter choices in general.

2.1 The (1+$\lambda$) EA

The first to conduct a rigorous runtime analysis of the (1+$\lambda$) EA were Jansen, De Jong, and Wegener [30]. They proved, among other results, that when optimizing OneMax a linear speed-up exists up to a population size of $\Theta(\log(n)\log\log(n)/\log\log\log(n))$, that is, for $\lambda=O(\log(n)\log\log(n)/\log\log\log(n))$, finding the optimal solution takes an expected number of $\Theta(n\log(n)/\lambda)$ generations, whereas for larger $\lambda$ at least $\omega(n\log(n)/\lambda)$ generations are necessary. This picture was completed in [13] with a proof that the expected number of generations taken to find the optimum is $\Theta(\frac{n\log n}{\lambda}+\frac{n\log\log\lambda}{\log\lambda})$. The implicit constants were determined in [27], giving the bound of $(1\pm o(1))(\frac{1}{2}\frac{n\ln\ln\lambda}{\ln\lambda}+\frac{e^{r}}{r}\frac{n\ln n}{\lambda})$, for any constant $r$, as mentioned in the introduction.

Aside from the optimization behavior on OneMax, not too much is known for the (1+$\lambda$) EA, or is at least not made explicit (it is easy to see that waiting times for an improvement which are larger than $\lambda$ reduce by a factor of $\Theta(\lambda)$ compared to one-individual offspring populations). Results made explicit are the $\Theta(n^{2}/\lambda+n)$ expected runtime (number of generations) on LeadingOnes [30], the worst-case $\Theta(n+n\log(n)/\lambda)$ expected runtime on linear functions [13], and the $O(m^{2}(\log n+\log w_{\max})/\lambda)$ runtime estimate for minimum spanning trees valid for $\lambda\leq m^{2}/n$ [39], where $n$ denotes the number of vertices of the input graph, $m$ the number of edges, and $w_{\max}$ the maximum of the integral edge weights.

2.2 Dynamic Parameter Choices

While it is clear the EAs with parameters changing during the run of the algorithm (dynamic parameter settings) can be more powerful than those only using static parameter settings, only recently considerable advantages of dynamic choices could be demonstrated by mathematical means (for discrete optimization problems; in continuous optimization, step size adaptation is obviously necessary to approach arbitrarily closely a target point). To describe the different ways to dynamically control parameters, we use in the following the language proposed in Eiben, Hinterding, and Michalewicz [26] and its extension from [12].

2.2.1 Deterministic Parameter Control

In this language, deterministic parameter control means that the dynamic choice of a parameter does not depend on the fitness landscape. The first to rigorously analyze a deterministic parameter control scheme are Jansen and Wegener [29]. They regard the performance of the (1+1) EA which uses in iteration $t$ the mutation rate $k/n$, where $k\in\{1,2,\ldots,2^{\lceil\log_{2}n\rceil-2}\}$ is chosen such that $\log_{2}(k)\equiv t-1\,\,\,(\text{mod}\,\lceil\log_{2}n\rceil-1)$. In other words, they cyclically use the mutation rates $1/n,2/n,4/n,\ldots,K/n$, where $K$ is the largest power of two that is less than $n/2$. Jansen and Wegener demonstrate that there exists an example function where this dynamic EA significantly outperforms the (1+1) EA with any static mutation rate. However, they also observe that for many classic problems, this EA is slower by a factor of $\Theta(\log n)$.

In [40], a rank-based mutation rate was analyzed for the ($\mu$+1) EA. A previous experimental study [6] suggested that this is a profitable approach, but the mathematical runtime analysis in [40] rather indicates the opposite. While there are artificial examples where a huge runtime gain could be shown and also the worst-case runtime of the ($\mu$+1) EA reduces from essentially $n^{n}$ to $O(3^{n})$, a rigorous analysis on the OneMax function rather suggests that the high rate of offspring generated with a mutation rate much higher than $1/n$ brings a significant risk of slowing down the optimization process.

For two non-standard settings in evolutionary computation, deterministic parameter control mechanisms also gave interesting results. For problems where the solution length is not known [5], more precisely, where the number or the set of bits relevant for the solution quality is unknown, again random mutation rates gave good results [16, 21]. Here however, not a power-law scheme was used, but rather one based on very slowly decreasing summable sequences. For problems where the discrete variables take many values, e.g., the search space is $\{0,\dots,r-1\}^{n}$ for some large $r$, the question is how to change the value of an individual variable. The results in [17] suggest that a harmonic mutation strength, that is, changing a variable value by $\pm i$ with $i$ chosen randomly with probability proportional to $1/i$, can be beneficial. This distribution was analyzed earlier in [8] for the one-dimensional case, where it was also shown to give the asymptotically best performance on a OneMax type problem.

For randomized search heuristics outside evolutionary computation, Wegener [43] showed that simulated annealing (using a time-dependent temperature) can beat the Metropolis algorithm (using a static temperature).

2.2.2 Adaptive Parameter Control

A parameter control scheme is called adaptive if it used some kind of feedback from the optimization process. This can be functionally dependent (e.g., the mutation rate depends on the fitness of the parent) or success-based (e.g., a 1/5th rule).

The first to conduct a runtime analysis for an adaptive parameter control mechanism (and show a small advantage over static choices) were Böttcher, Doerr, and Neumann [3]. They proposed to use the fitness-dependent mutation rate of $1/(\textsc{LeadingOnes}(x)+1)$ for the optimization of the LeadingOnes test function. They proved that with this choice, the runtime of the (1+1) EA improves to roughly $0.68n^{2}$ compared to a time of $0.86n^{2}$ stemming from the classic mutation rate $1/n$ or a runtime of $0.77n^{2}$ stemming from the asymptotically optimal static rate of approximately $1.59/n$.

For the $(1+(\lambda,\lambda))$ GA, a fitness-dependent offspring population size of order $\lambda=\Theta(\sqrt{n/d(x)}\,)$ was suggested in [15], where $d(x)$ is the fitness-distance of the parent individual to the optimum. This choice improves the optimization time (number of fitness evaluations until the optimum is found) on OneMax from $\Theta(n\sqrt{\log(n)\log\log\log(n)/\log\log(n)}\,)$ stemming from the optimal static parameter choice [10] to $O(n)$. Since in this adaptive algorithm the mutation rate $p$ is functionally dependent on the offspring population size, namely via $p=\lambda/n$, the dynamic choice of $\lambda$ is equivalent to a fitness-dependent mutation rate of $1/\sqrt{nd(x)}$.

In the aforementioned work by Badkobeh et al. [2], a fitness-dependent mutation rate of \max\bigl{\{}\frac{\ln\lambda}{n\ln(en/d(x))},\frac{1}{n}\bigr{\}} was shown to improve the classic runtime of O\bigl{(}\frac{n\log\log\lambda}{\log\lambda}+\frac{n\log n}{\lambda}\bigr{)} to O\bigl{(}\frac{n}{\log\lambda}+\frac{n\log n}{\lambda}\bigr{)}. In [19], the (1+1) EA using a $k$-bit flip mutation operator together with a fitness-dependent choice of $k$ was shown to give a performance on OneMax that is very close to the theoretical optimum (among all unary unbiased black-box algorithms), however, this differs only by lower order terms from the performance of the simple randomized local search heuristic (RLS). For nature-inspired algorithms other than evolutionary ones, Zarges [44, 45] proved that fitness-dependent mutation rates can be beneficial in artificial immune systems.

2.2.3 Self-adjusting and Self-adaptive Parameter Control

While all these results show an advantage of adaptive parameter settings, it remains questionable if an algorithm user would be able to find such a functional dependence of the parameter on the fitness. This difficulty can be overcome via self-adjusting parameter choices, where the parameter is modified according to a simple rule often based on the success or progress of previous iterations, or via self-adaptation, where the parameter is encoded in the genome and thus subject to variation and selection. The understanding of self-adaptation is still very limited. The only theoretical work on this topic [7], however, is promising and shows examples where self-adaptation can lead to significant speed-ups for non-elitist evolutionary algorithms.222Note added in proof: Very recently, it was shown in [25] that the (1,$\lambda$) EA with self-adaptive mutation rate has a runtime asymptotically equivalent to the one of the algorithm regarded in this work.

In contrast to this, the last years have produced a profound understanding of self-adjusting parameter choices. The first to perform a mathematical analysis were Lässig and Sudholt [34], who considered the (1+$\lambda$) EA and a simple parallel island model together with two self-adjusting mechanisms for population size or island number, including halving or doubling it depending on whether the current iteration led to an improvement or not. These mechanisms were proven to give significant improvements of the “parallel” runtime (number of generations) on various test functions without increasing significantly the “sequential” runtime (number of fitness evaluations).

In [12] it was shown that the fitness-dependent choice of $\lambda$ for the $(1+(\lambda,\lambda))$ GA described above can also be found in a self-adjusting way. To this aim, another success-based mechanism was proposed, which imitates the $1/5$-th rule from evolution strategies. With some modifications, this mechanism also works on random satisfiability problems [4]. For the problem of optimizing an $r$-valued OneMax function, a self-adjustment of the step size inspired by the $1/5$-th rule was found to find the asymptotically best possible runtime in [18].

These results indicate that success-based dynamics work well for adjusting parameters when a monotonic relation like “if progress is difficult, then increase the population size” holds. For adjusting a parameter like the mutation rate, it is less obvious how to do this. For example, in the search space $\{0,1\}^{n}$ both a too large mutation rate (creating a stronger drift towards a Hamming distance of $n/2$ from the optimum) and a too small mutation rate (giving a too small radius of exploration) can be detrimental. For this reason, to obtain a self-adjusting version of their result on the optimal number $k$ to optimize OneMax via $k$-bit flips [19], in [20] a learning mechanism was proposed that from the medium-term past estimates the efficiency of different parameter values. As shown there, this does find the optimal mutation strength sufficiently well to obtain essentially the runtime stemming from the fitness-dependent mutation strength exhibited before.

In the light of these works, our result from the methodological perspective shows that some of the difficulties of the learning mechanism of [20], e.g., the whole book-keeping being part of it and also the setting of the parameters regulating how to discount information over time, can be overcome by the mechanism proposed in this work. In a sense, the use of larger populations enables us to adjust the mutation rate solely on information learned in the current iteration. However, we do also use the idea of [20] to intentionally use parameter settings which appear to be slightly off the current optimum to gain additional insight.

2.2.4 Selection Hyper-Heuristics

We note that so-called selection hyper-heuristics may lead to processes resembling dynamic parameter choices. Selection hyper-heuristics are methods that select, during the run of the algorithm, which one out of several pre-specified simpler algorithmic building blocks to use. When the different pre-specified choices are essentially identical apart from an internal parameter, then this selection hyper-heuristic approach could also be interpreted as a dynamic choice of the internal parameter. For example, when only the two mutation operators are available that flip exactly one or exactly two bits, then a selection hyper-heuristic choosing between them could also be interpreted as the randomized local search heuristic using a dynamic choice of the number of bits it flips. We briefly review the main runtime analysis works of this flavor.

The first to conduct a rigorous runtime analysis for selection hyper-heuristics were Lehre and Öczan [35]. They show that the (1+1) EA using the mixed strategy of choosing the mutation operator randomly between the 1-bit flip operator (with probability $p$) and the $2$-bit flip operator (with probability $1-p$) optimizes the OneMax function in time $\min\{\frac{1}{p}\,n(1+\ln n),\,\frac{1}{1-p}\,n^{2}(1-\frac{1}{n})\}$. It appears to us that, most likely, this result is not absolutely correct, since, e.g., in the case $p=0$ the expected optimization time is clearly infinite (if the random initial search point has an odd Hamming distance from the optimum, then the optimum cannot be reached only via $2$-bit flips). So most likely, the upper bound should be increased by an additive term of $\frac{1}{p}n$.

Lehre and Öczan [35] further construct an example function GapPath, which has the property that it can only be optimized when both the $1$-bit and the $2$-bit flip mutation operator occur with positive probability and use this to discuss various static and dynamic ways to randomly decide between the two operators. A similar example for a multi-objective optimization setting was given by Qian, Tang, and Zhou [41].

Note that for $p=\frac{1}{2}$, the randomized selection heuristic is the classic randomized search heuristic using $1$-bit and $2$-bit flips, which was regarded, among others, in [28, 39].

In [1], Alanazi and Lehre extent the previous work to several classic selection hyper-heuristics like simple random (take a random low-level heuristic in each iteration), random gradient (take a random low-level heuristic and repeat using it as long as a true fitness improvement is obtained), greedy (in each iteration, use all low-level heuristic in parallel and continue with best result obtained), permutation (generate initially a random cyclic order of the low-level heuristics and then use them in this order). Again for the choice between $1$-bit and $2$-bit flips, they prove upper and lower bounds for the expected optimization time on the LeadingOnes benchmark function. While the results are relatively tight (the corresponding upper and lower bounds deviate by at most a factor of $6$), the intervals of possible asymptotic runtimes intersect, so this first runtime analysis of this type does not yet give a conclusive picture on how the different selection heuristics compare.

Given that the probabilities to find a true improvement are very low in this discrete optimization problem, one would expect that all these four heuristics have the same runtime (apart from lower order terms), and this is indeed the first set of results by Lissovoi, Oliveto, and Warwicker [37], who show that the expected runtime in all cases is $(1+o(1))\frac{1}{2}\ln(3)n^{2}$. They build on this strong result by proposing to use a slower adaptation. For the random gradient method, they propose to use a random low-level heuristic for up to $\tau$ iterations. If an improvement is found, immediately another phase with this operator is started. If a phase of $\tau$ iterations does not see a fitness improvement, then a new phase is started with a random operator.

For this generalized random gradient heuristic using a phase length of $\tau=cn$ for a constant $c$, they show (still for the LeadingOnes problem and the $1$-bit and $2$-bit mutation operator) an expected runtime of $(1+o(1))g(c)n^{2}$, where $g(c)$ is a constant depending on $c$ only that tends to $\frac{\ln(2)+1}{4}$ when $c$ is tending to infinity. Consequently, this new hyper-heuristic outperforms the previously investigated ones and approaches arbitrarily well the best-possible performance that can be obtained from the two mutation operators, which is, as also shown in [37], $(1+o(1))\frac{\ln(2)+1}{4}n^{2}$.

The generalized random gradient heuristic was further extended in [24]. There an operator was defined as successful (which leads to another phase using this operator) if it leads to $\sigma$ improvements in a phase of at most $\tau$ iterations. Hence in this language, the previous generalized random gradient heuristic uses $\sigma=1$. By using a larger value of $\sigma$, the algorithm is able to take more robust decisions on what is a success. This was used in [24] to determine the phase length $\tau$ in a self-adjusting manner. While the previous work [37] does not state this explicitly, the choice of $\tau$ is crucial. A $\tau$-value of smaller asymptotic order than $\Theta(n)$ leads to typically no improvement within a phase and thus reverts the algorithms to the simple random heuristic. A $\tau$-value of more than $cn\ln(n)$, where $c$ is a suitable constant, results in that both operators are successful in most parts of the search space. Consequently, the algorithm sticks to the first choice for a large majority of the optimization process and thus does not profit from the availability of both operators.

Since the choice of $\tau$ is that critical, a mechanism successfully adjusting it to the right value is desirable. In [24] it is shown that by choosing $\sigma\in\Omega(\log^{4}n)\cap o(\sqrt{n/\log n})$—note that this is a quite wide range—the value of $\tau$ can be easily adjusted on the fly via a multiplicative update rule resembling vaguely the $1/5$-th success rule. This gives again the asymptotically optimal runtime of $(1+o(1))\frac{\ln(2)+1}{4}n^{2}$.

3 Preliminaries

We shall now formally define the algorithm analyzed and present some fundamental tools for the analysis.

3.1 Algorithm

We consider the (1+$\lambda$) EA with self-adjusting mutation rate for the minimization of pseudo-boolean functions $f:\{0,1\}^{n}\to\mathbb{R}$, defined as Algorithm 1.

The general idea of the mutation scheme is to adjust the mutation strength according to its success in the population. We perform mutation by applying standard bit mutation with two different mutation probabilities $r/(2n)$ and $2r/n$ and we call $r$ the mutation rate. More precisely, for an even number $\lambda\geq 2$ the algorithm creates $\lambda/2$ offspring with mutation rate $r/2$ and with $2r$ each.

The mutation rate is adjusted after each selection. With probability a half, the new rate is taken as the mutation rate that the best individual (i. e. the one with the lowest fitness, ties broken uniformly at random) was created with (success-based adjustment). With the other 50% probability, the mutation rate is adjusted to a random value in $\{r/2,2r\}$ (random adjustment). Note that the mutation rate is adjusted in each iteration, that is, also when all offspring are worse than the parent and thus the parent is kept for the next iteration.

If an adjustment of the rate results in a new rate $r$ outside the interval $[2,n/4]$, we replace this rate with the corresponding boundary value. Note that in the case of $r<2$, a subpopulation with rate less than $1$ would be generated, which means flipping less than one bit in expectation. At a rate $r>n/4$, a subpopulation with rate larger than $n/2$ would be created, which again is not a very useful choice.

We formulate the algorithm to start with an initial mutation rate $r^{\operatorname{init}}$. The only assumption on $r^{\operatorname{init}}$ is to be greater than or equal to $2$. The (1+$\lambda$) EA with this self-adjusting choice of the mutation rate is given as pseudocode in Algorithm 1.

Let us explain the motivation for the random adjustments of the rate. Without such random adjustments, the rate can only be changed into some direction if a winning offspring is generated with this rate. For simple functions like OneMax, this is most likely sufficient. However, when the fitness of the best of $\lambda/2$ offspring, viewed as a function of the rate, is not unimodal, then several adjustments into a direction at first not yielding good offspring might be needed to reach good values of the rate. Here, our random adjustments enable the algorithm to cross such a valley of unfavorable rate values. We note that such ideas are not uncommon in evolutionary computation, with standard-bit mutation being the most prominent example (allowing to perform several local-search steps in one iteration to cross fitness valleys).

A different way to implement a mechanism allowing larger changes of the rate to cross unfavorable regions would have been to not only generate offspring with rates $r/2$ and $2r$, but to allow larger deviations from the current rate with some small probability. One idea could be choosing for each offspring independently the rate $r2^{-i}$ with probability $2^{-|i|-1}$ for all $i\in\mathbb{Z}$, $i\neq 0$. This should give similar results, but to us the process appears more chaotic (e.g., with not the same number of individuals produced with rates $r/2$ and $2r$).

The runtime, also called the optimization time, of the (1+$\lambda$) EA is the smallest $t$ such that an individual of minimum $f$-value has been found. Note that $t$ corresponds to a number of iterations (also called generations), where each generation creates $\lambda$ offspring. Since each of these offspring has to be evaluated, the number of function evaluations, which is a classical cost measure, is by a factor of $\lambda$ larger than the runtime as defined here. However, assuming a massively parallel architecture that allows for parallel evaluation of the offspring, counting the number of generations seems also a valid cost measure. In particular, a speed-up on the function $\textsc{OneMax}(x_{1},\dots,x_{n}):=x_{1}+\dots+x_{n}$ by increasing $\lambda$ can only be observed in terms of the number of generations. Note that for reasons of symmetry, it makes no difference whether OneMax is minimized (as in the present paper) or maximized (as in several previous research papers).

Throughout the paper, all asymptotic notation will be with respect to the problem size $n$.

3.2 Drift Theorems

Our results are obtained by drift analysis, which is also used in previous analyses of the (1+$\lambda$) EA without self-adaptation on OneMax and other linear functions [13, 27].

The first theorems stating upper bounds on the hitting time using variable drift go back to [31, 38]. We take a formulation from [36] but simplify it to Markov processes for notational convenience.

Theorem 3** (Variable Drift, Upper Bound).**

Let $(X_{t})_{t\geq 0}$, be random variables describing a Markov process over a finite state space $S\subseteq\{0\}\cup[x_{\mathrm{min}},x_{\mathrm{max}}]$, where $x_{\mathrm{min}}>0$. Let $T$ be the random variable that denotes the earliest point in time $t\geq 0$ such that $X_{t}=0$. If there exists a monotone increasing function $h(x)\colon[x_{\mathrm{min}},x_{\mathrm{max}}]\to\mathbb{R}^{+}$, where $1/h(x)$ is integrable on $[x_{\mathrm{min}},x_{\mathrm{max}}]$, such that for all $x\in S$ with $\operatorname{Pr}(X_{t}=x)>0$ we have

[TABLE]

then for all $x^{\prime}\in S$ with $\operatorname{Pr}(X_{0}=x^{\prime})>0$

[TABLE]

The variable drift theorem is often applied in the special case of additive drift in discrete spaces: assuming $E(X_{t}-X_{t+1}\mid X_{t}=x;X_{t}>0)\geq\varepsilon$ for some constant $\varepsilon$, one obtains $E(T\mid X_{0}=x^{\prime})\leq x^{\prime}/\varepsilon$.

Since we will make frequent use of it in the following sections as well, we will also give the version of the Multiplicative Drift Theorem for upper bounds, due to [14]. Again, this is implied by the previous variable drift theorem.

Theorem 4** (Multiplicative Drift [14]).**

Let $(X_{t})_{t\geq 0}$ be random variables describing a Markov process over a finite state space $S\subseteq\mathbb{R}^{+}_{0}$ and let $x_{\mathrm{min}}:=\min\{x\in S\;|\;x>0\}$. Let $T$ be the random variable that denotes the earliest point in time $t\geq 0$ such that $X_{t}=0$. If there exist $\delta>0$ such that for all $x\in S$ with $\operatorname{Pr}(X_{t}=x)>0$ we have

[TABLE]

then for all $x^{\prime}\in S$ with $\operatorname{Pr}(X_{0}=x^{\prime})>0$,

[TABLE]

3.3 Chernoff Bounds

For reasons of self-containedness, we state two well-known multiplicative Chernoff bounds and a lesser known additive Chernoff bound that is also known in the literature as Bernstein’s inequality. These bounds can be found, e.g., in [9].

Theorem 5** (Bernstein’s inequality, Chernoff bounds).**

Let $X_{1},\ldots,X_{n}$ be independent random variables and $X=\sum_{i=1}^{n}X_{i}$.

1. (a)

Let $b$ be such that $E(X_{i})-b\leq X_{i}\leq E(X_{i})+b$ for all $i=1,\ldots,n$. Let $\sigma^{2}=\sum_{i=1}^{n}\mathrm{Var}(X_{i})=\mathrm{Var}[X]$. Then for all $\Delta\geq 0$,

[TABLE] 2. (b)

Assume that for all $i=1,\dots,n$, the random variable $X_{i}$ takes values in $[0,1]$ only. Then

• •

$\operatorname{Pr}(X\leq(1-\delta)E(X))\leq\exp(-\delta^{2}E(X)/2)$* for all $0\leq\delta\leq 1$,*

• •

$\operatorname{Pr}(X\geq(1+\delta)E(X))\leq\exp(-\delta^{2}E(X)/(2+\delta))$* for all $\delta>0$.*

3.4 Occupation Probabilities

As mentioned above, we will be analyzing two depending stochastic processes: the random decrease of fitness and the random change of the mutation rate. Often, we will prove by drift analysis that the rate is drifting towards values that yield an almost-optimal fitness decrease. However, once the rate has drifted towards such values, we would also like the rates to stay in the vicinity of these values in subsequent steps. To this end, we apply the following theorem from [33]. Note that in the paper a slightly more general version including a self-loop probability is stated, which we do not need here.

Theorem 6** (Theorem 7 in [33]).**

Let a Markov process $(X_{t})_{t\geq 0}$ on $\mathbb{R}_{0}^{+}$, where $\lvert X_{t}-X_{t+1}\rvert\leq c$, with additive drift of at least $d$ towards [math] be given (i. e., $E(X_{t}-X_{t+1}\mid X_{t};X_{t}>0)\geq d$), starting at [math] (i.e. $X_{0}=0$). Then we have, for all $t\in\mathbb{N}$ and $b\in\mathbb{R}_{0}^{+}$,

[TABLE]

We can readily apply this theorem in the following lemma that will be used throughout the paper to bound the rate $r_{t}$.

Lemma 7**.**

If there is a point $a\geq 4$ such that $\operatorname{Pr}(r_{t+1}<r_{t}\mid r_{t}>a)\geq 1/2+\varepsilon$ for some constant $\varepsilon>0$, then for all $t^{\prime}\geq\min\{t\mid r_{t}\leq a\}$ and all $b\geq 2$ it holds $\operatorname{Pr}(r_{t^{\prime}}\geq a\cdot 2^{b+1})\leq 2e^{-2b\varepsilon/3}$.

Proof.

Apply Theorem 6 to the process $X_{t}:=\max\{0,\lfloor\log_{2}(r_{t}/a)\rfloor\}$. Note that this process is on $\mathbb{N}_{0}$, moves by an absolute value of at most $1$ and has drift $E(X_{t}-X_{t+1}\mid X_{t};X_{t}>0)\geq 2\varepsilon$. We use $c:=1$ and $d:=2\varepsilon$ in the theorem and estimate $1-b\leq-b/2$. ∎

3.5 Useful Inequalities

The following well-known estimates will be used regularly in this work.

Lemma 8**.**

1. (a)

For all $x>0$, we have $\ln x\leq x/e$. 2. (b)

For all $x\in\mathbb{R}$, we have $1+x\leq e^{x}$. 3. (c)

For all $0\leq x\leq 2/3$, we have $1-x\geq e^{-x-x^{2}}\geq e^{-2x}$.

Proof.

Since $(\ln x-x/e)^{\prime}=1/x-1/e,x>0$, then $(e,0)$ is the maximum point of $\ln x-x/e,x>0$. Thus part (a) holds.

We notice that $(1+x-e^{x})^{\prime}=1-e^{x}$, then $(0,0)$ is the maximal point of $1+x-e^{x}$. Therefore part (b) holds.

For the proof of part (c), the second inequality can be easily obtained, since $e^{x}$ monotonically increases and $-x-x^{2}\geq-2x$ when $0\leq x\leq 2/3$. To prove the first inequality, let $f(x)=1-x-e^{-x-x^{2}}$ for $0\leq x\leq 2/3$. Then $f^{\prime}(x)=(1+2x)e^{-x-x^{2}}-1$ and $f^{\prime\prime}(x)=(1-4x-4x^{2})e^{-x-x^{2}}$. We notice that $f^{\prime\prime}(x)$ has only one zero point between [math] and $2/3$. Let $x_{0}$ denote the zero point. Thus $f^{\prime\prime}(x)$ is positive before $x_{0}$, and then becomes negative after $x_{0}$. Using the fact $f^{\prime}(0)=0$ and $f^{\prime}(2/3)<0$, we can easily see that $f^{\prime}(x)$ first increases from [math] to $f^{\prime}(x_{0})>0$, and then decreases below [math]. This means the minimum of $f(x)$ is attained at [math] or $2/3$. Due to the fact that $f(0)=0$ and $f(2/3)>0$, the first inequality in part (c) holds. ∎

4 Proof Overview

Since the following runtime analysis of our self-adjusting (1+$\lambda$) EA on the OneMax function is slightly technical, let us outline the main proof ideas here in an informal manner. Let always $x$ denote the parent individual of the current iterations and $k:=f(x)$ its fitness distance from the optimum (recall that we are minimizing the OneMax function, hence the fitness distance equals the objective function value which in turn is the Hamming distance from the optimum $x^{*}=(0,\dots,0)\in\{0,1\}^{n}$ and thus the number of ones in $x$).

The outer proof argument is variable drift, that is, for each fitness value $k$ we estimate the expected fitness gain (“progress”) in an iteration starting with a parent individual $x$ with $f(x)=k$ and then we use the variable drift theorem (Theorem 3) to translate this information on the progress into an expected runtime (number of generations until the optimum is found).

To obtain sufficiently strong lower bounds on the expected progress, we need to argue that the rate self-adjustment sufficiently often sets the current rate to a sufficiently good value. This is the main technical difficulty as it needs a very precise analysis of the quality of the offspring in both subpopulations. Since this requires different arguments depending on the current fitness distance $k$, we partition the process into three regimes.

In the far region covering the fitness distance values $k\geq n/\ln(\lambda)$, we need a rate $r$ that is at least almost logarithmic in $\lambda$ to ensure that the average progress is high enough. Note that to gain the $\Theta(n)$ fitness levels from the initial value of approximately $n/2$ to $n/\ln\lambda$ in at most a time of $O(n/\log\lambda)$, we need an average progress of $\Omega(\ln\lambda)$ per iteration. We shall not be able to obtain this progress throughout this region, but we will obtain an expected progress of

[TABLE]

in an iteration with initial fitness distance $k$. This will be sufficient to reach a fitness distance of $n/\log\lambda$ in the desired $O(n/\log\lambda)$ iterations.

As said, the main difficulty is arguing that the self-adjusting mechanism keeps the rate sufficiently often in the range we need, which is roughly

[TABLE]

for all $k\in[n/\ln\lambda,n/2-1]$. Note that these range boundaries, in particular the upper one, depend strongly on $k$. Hence for arguing that our rate is sufficiently often in this range, we need to consider both the rate changes from the self-adjustments and the changing $k$-value. In this region we profit from the fact that we work with relatively high rates, which lead to strong concentration behaviors. This allows to argue that, for example, exceeding the upper boundary already by small constant factors is highly unlikely.

The middle region covering the fitness distances $k\in[n/\lambda,n/\ln\lambda]$ is small enough so that we do not require the algorithm to find a near-optimal rate very often. In fact, it suffices that the rate is below $\frac{1}{2}\ln\lambda$ with constant probability. This ensures an expected progress of at least $\min\{\frac{1}{8},\sqrt{\lambda}k/32n\}$, which is sufficient to traverse also this region in time $O(n/\log\lambda)$. Consequently, in this region we only need to argue that the rate does not become too large, whereas there is no lower bound on $r$ which we require. Further, the upper bound of $\tfrac{1}{2}\ln\lambda$ is large enough so that strong concentration arguments can be exploited, which show that deviations above the upper bound are highly unlikely. Since the upper bound does not change over time, we can now conveniently use an occupation probability argument to show that the rate with constant probability is only a small constant factor over the desired range, which is enough since the random fluctuations of the rate with constant probability reduce the rate further to the desired range.

In the near region covering the remaining fitness distance values $k\leq n/\lambda$, the parent individual is already so close to the optimum that any rate higher than the minimal rate $r=2$ is sub-optimal. Hence in this region, we know precisely the optimum value for the rate. Nevertheless, it is not very easy to show that the subpopulation using the smaller rate $r/2$ has a higher chance of containing the new parent individual. Due to the small rates predominant in this region and the fact that progress generally is difficult (due to the proximity to the optimum), often both subpopulations will contain best offspring (which are in fact copies of the parent). Hence the typical reason for the winning offspring stemming from the smaller-rate subpopulation is not anymore that the other subpopulation contains only worse individuals, but only that the smaller-rate subpopulation contains more best offspring.

By quantifying this effect we establish a drift of the rate down to its minimum value $2$ in the near region, and another occupation probability argument is used to show that the rate with sufficiently high probability will take this value in subsequent sets. Nevertheless, the concentration of the rate around this target value is weaker than in the other regions, and additional care has to be taken to handle iterations in the near region in which the unlikely, but still possible event occurs that the rate exceeds $\ln\lambda$. This value is a critical threshold in the near region since rates larger than $\ln\lambda$ will typically make all offspring worse than the parent. As an additional obstacle, there is a small interval of rates above the threshold in which we cannot show a drift of the rate to its minimum. To show that this small interval nevertheless is left towards its lower end with probability $\Omega(1)$, the occupation probability argument is combined with a potential function whose shape exploits the random adjustments that the (1+$\lambda$) EA performs with 50% probability.

5 Far Region

In this first of three technical sections, we analyze the optimization behavior of our self-adjusting (1+$\lambda$) EA in the regime where the fitness distance $k$ is at least $n/\ln\lambda$. Since we are relatively far from the optimum, it is relatively easy to make progress. On the other hand, this regime spans the largest number of fitness levels (namely $\Theta(n)$), so we need to exhibit a sufficient progress in each iteration. Also, this is the regime where the optimal mutation rate varies most drastically. Since it is not important for the following analysis, we remark without proof that the optimal rate is $n$ for $k\geq n/2+\omega(\sqrt{n\log\lambda})$, $(1+o(1))n/2$ for $k=n/2\pm o(\sqrt{n\log\lambda})$, and then quickly drops to $r=\Theta(\log\lambda)$ for $k\leq n/2-\varepsilon n$. Despite these difficulties, our (1+$\lambda$) EA manages to find sufficiently good mutation rates to be able to reach a fitness distance of $k=n/\ln\lambda$ in an expected number of $O(n/\log\lambda)$ iterations.

Lemma 9**.**

Let $n$ be sufficiently large and $0<k<n/2$. We define $c_{1}(k)=(2\ln(en/k))^{-1}$ and $c_{2}(k)=4n^{2}/(n-2k)^{2}$.

1. (a)

If $n/\ln\lambda\leq k$ and $r\leq c_{1}(k)\ln\lambda$, then the probability that a best offspring has been created with rate $2r$ is at least $0.64$. 2. (b)

Let $\lambda\geq 100$. If $c_{2}(k)\ln\lambda\leq r\leq n/4$, then the probability that all best offspring have been created with rate $r/2$ is at least $0.51$. 3. (c)

If $r\geq 2(1+\gamma)c_{2}(k)\ln\lambda$, then the probability that all best offspring are worse than the parent is at least $1-\lambda^{-\gamma}$.

Proof.

To prove part (a), let $q(k,i,r)$ and $Q(k,i,r)$ be the probabilities that standard bit mutation with mutation rate $p=r/n$ creates from a parent with fitness distance $k$ an offspring with fitness distance exactly $k-i$ and at most $k-i$. Then

[TABLE]

and $Q(k,i,r)=\sum_{j=i}^{k}q(k,j,r)$. We first show that the probability of not achieving $i\geq 2r$ is less than $0.2$. This is because for large enough $n$ and $r\geq 2$ we have $(1-o(1))\binom{k}{2r}\geq(1-o(1))4!(k/(2r))^{2r}>4(k/(2r))^{2r}$ and for $r\leq c_{1}(k)\ln\lambda=o(\sqrt{n})$ by Lemma 8, part (c) we have $(1-2r/n)^{n}\geq e^{-2r-4r^{2}/n}=(1-o(1))e^{-2r}$, thus

[TABLE]

Considering $\lambda/2$ offspring using rate $2r$, the probability that none of them achieves a fitness improvement of at least $2r$ is less than $(1-4/\lambda)^{\lambda/2}<\exp(-4/2)<0.2$. We next argue that an offspring having a progress of $2r$ or more with good probability comes from the $2r$-subpopulation. We notice that

[TABLE]

Thus for $r\geq 2$ and $i\geq 2r$ we have

[TABLE]

Therefore if the best progress among all $\lambda$ offspring is $i$ and $i\geq 2r$, the conditional probability that an offspring having fitness distance $k-i$ is generated with rate $2r$ is at least $q(k,i,2r)/(q(k,i,2r)+q(k,i,r/2))\geq 4/5$. In total, the probability that a best offspring has been created with rate $2r$ is at least

[TABLE]

For the proof of part (b), let $\lambda\geq 100$ and $c_{2}(k)\ln\lambda\leq r\leq n/4$. Our idea is to show a probability of $o(1/\lambda)$ for the event that an offspring with rate $2r$ is not worse than the expected fitness of an offspring with rate $r/2$. Let $X(k,r)$ denote the random decrease of the fitness distance when applying standard bit mutation with probability $p=r/n$ to an individual with $k$ ones. Then

[TABLE]

According to Bernstein’s inequality (Theorem 5 (a)), for any $\Delta>0$ we have

[TABLE]

We apply this bound with $X=X(k,2r)$ and $\Delta=E(X(k,r/2))-E(X(k,2r))=(n-2k)1.5r/n>0$. Then,

[TABLE]

We notice that we have $7n-8r-6k>7n-4n-3n=0$ in the second inequality. From a union bound we see that with probability less than $\lambda^{-9/7}(\lambda/2)<100^{-2/7}/2<0.14$, the best offspring created with rate $2r$ is at least as good as the expected fitness of an individual created with rate $r/2$.

We estimate the probability that the best offspring using rate $r/2$ has a fitness distance at most $E(X(k,r/2))$. Let $y$ be an offspring obtained from $x$ via standard bit mutation with mutation rate $r/(2n)$. Let $X^{+}$ and $X^{-}$ be the number of one-bits flipped and zero-bits flipped, respectively.

[TABLE]

Both $X^{+}$ and $X^{-}$ are binomially distributed and $X(k,r/2)=X^{+}-X^{-}$. We aim at a lower bound for $\Pr(X(k,r/2)\geq E(X(k,r/2)))$. We have $\Pr(X^{+}\geq E(X^{+})-1)\geq 1/2$, since the median of $X^{+}$ is between $\lfloor E(X^{+})\rfloor$ and $\lceil E(X^{+})\rceil$ by [32]. It remains to bound $\operatorname{Pr}(X^{-}\leq E(X^{-})-1)$. We notice that $E(X^{-})=(n-k)r/(2n)>\ln\lambda>1$ and $E(X^{-})\leq(n-k)/4$. Applying Theorem 12 in [11] to binomial random variable $\tilde{X}^{-}:=(n-k)-X^{-}$, we obtain $\operatorname{Pr}(X^{-}\leq E(X^{-})-1)=\operatorname{Pr}(\tilde{X}^{-}\geq E(\tilde{X}^{-})+1)\geq 0.037.$ Therefore

[TABLE]

For $\lambda/2$ offspring using rate $r/2$, the probability that the best one has a fitness distance at most $E(X(k,r/2))$ is more than $1-(1-0.0185)^{\lambda/2}\geq 1-0.9815^{50}\geq 0.6$. Therefore with probability at least $0.6\cdot(1-0.14)>0.51$, all best offspring are from the $r/2$-subpopulation. This proves the second statement of the lemma.

For proof of part (c), let $r\geq 2(1+\gamma)c_{2}(k)\ln\lambda$. An offspring created with mutation rate $r/n$ is at least as good as its parent if and only if $X(k,r)\geq 0$. By using Bernstein’s inequality (Theorem 5 (a)) with $X=X(k,r)$ and $\Delta=-E(X(k,r))$, we have

[TABLE]

Since $r\geq 2(1+\gamma)c_{2}(k)\ln\lambda$, the corresponding probabilities for rate $r/2$ and $2r$ are at most $1/\lambda^{1+\gamma}$ and $1/\lambda^{4+4\gamma}$, respectively. By a union bound, with probability at most $1/\lambda^{\gamma}$, the best offspring is at least as good as its parent. This proves part (c). ∎

The lemma above shows that the rate $r$ is subject to a constant drift towards the interval $[c_{1}(k)\ln n,c_{2}(k)\ln n]$. Unfortunately, we cannot show that we obtain a sufficient fitness progress for all $r$-values in this range. However, we can do so for a range smaller only by constant factors. This is what we do now (for large values of $k$) and in Lemma 11 (for smaller values of $k$). This case distinction is motivated by the fact that $c_{2}(k)$ becomes very large when $k$ approaches $n/2$. Having a good fitness drift only for such a smaller range of $r$-values is not a problem since the random movements of $r$ let us enter the smaller range with constant probability. This is what we will exploit in Theorem 12 and its proof.

Let $2\leq r\leq n/4$ be the current rate and let $\tilde{r}\in\{r/2,2r\}$. Let $\tilde{\Delta}(\lambda/2,k,\tilde{r})$ denote the fitness gain of the best of $\lambda/2$ offspring generated with rate $\tilde{r}$ from a parent $x$ with fitness distance $k:=\textsc{OneMax}(x)$ and the parent itself. More precisely, let $x^{(i)},i\in\{1,\dots,\lambda/2\}$, be independent offspring generated from $x$ by flipping each bit independently with probability $\tilde{r}/n$. Then the random variable $\tilde{\Delta}(\lambda/2,k,\tilde{r})$ is defined by $\max\{0,k-\min\{\textsc{OneMax}(x^{(i)})\mid i\in\{1,\dots,\lambda/2\}\}\}$. We use $\Delta:=\max\{\tilde{\Delta}(\lambda/2,k,r/2),\\ \tilde{\Delta}(\lambda/2,k,2r)\}$ to denote the fitness gain in a iteration which uses $x$ as parent and $r$ as mutation rate.

We next show that a region contained in $[c_{1}(k)\ln\lambda,c_{2}(k)\ln\lambda]$ provides at least a logarithmic (in $\lambda$) drift on fitness.

Lemma 10**.**

Let $n$ be sufficiently large, $2n/5\leq k<n/2$ and $\lambda\geq e^{5}>148$. If $2\leq\ln(\lambda)\leq r\leq\min\{n^{2}\ln(\lambda)/(25(n-2k)^{2}),n/4\}$, then $E(\Delta\mid k)\geq 10^{-3}\ln\lambda$.

Proof.

We first notice that $\ln(\lambda)\geq 5$ and $n^{2}/(25(n-2k)^{2})\geq 1/(25\cdot 0.2^{2})=1$ for all $2n/5\leq k<n/2$. We aim to prove $\tilde{\Delta}(\lambda/2,k,\tilde{r})\geq 10^{-3}\ln\lambda$ with $\tilde{r}=r/2$. Let $X^{+}$ and $X^{-}$ be the number of one-bits flipped and zero-bits flipped, respectively, in an offspring using rate $p=\tilde{r}/n$. $X^{+}$ and $X^{-}$ follow binomial distributions $\operatorname{Bin}(k,p)$ and $\operatorname{Bin}(n-k,p)$, respectively. Let $u:=E(X^{+})=kp$. Then $u\leq k/8$ and $u\geq(k/n)\ln(\lambda)/2\geq 0.4\cdot 0.5\cdot\ln\lambda=0.2\ln\lambda\geq 1$. Furthermore, let

[TABLE]

Applying Theorem 10 in [11], we obtain $F(u)=\Pr(X^{+}\geq E(X^{+}))\geq 1/4$. Similarly $\Pr(X^{-}\leq E(X^{-}))\geq 1/4$. We prove that for $\delta:=E(X^{-})-E(X^{+})+0.05\ln\lambda=(n-2k)p+0.05\ln\lambda$, we have $F(u+\delta)\geq\lambda^{-0.8}/9$, and thus $\Pr(X^{+}-X^{-}\geq 0.05\ln\lambda)\geq\lambda^{-0.8}/36$. Since for any $x\in\mathbb{Z}_{\geq u}$ we have

[TABLE]

we obtain $B(\lceil u\rceil)>B(\lceil u\rceil+1)>\dots>B(k)$, and thus $F(u+\delta)\geq\lceil\delta\rceil B(\lceil u+2\delta\rceil)$ as well as $\lceil\delta\rceil B(\lceil u\rceil)\geq F(u)-F(u+\delta)$. We see that

[TABLE]

where we used $k(1-p)=k-u\geq 8u-u=7u$. Using a sharp version of Stirling’s approximation due to Robbins [42], we compute

[TABLE]

We notice that $u\geq(k/n)\ln(\lambda)/2\geq 0.4\cdot 0.5\cdot\ln\lambda=0.2\ln\lambda$ and

[TABLE]

Thus $2\delta/(3u)\leq 2/3$ and $2\delta/(u+2\delta)\leq 2\delta/(\delta+2\delta)=2/3$. By Lemma 8 (c) we have $\ln(1-x)\geq-x-x^{2}\geq-2x$ for $0\leq x\leq 2/3$. Hence

[TABLE]

where in the last inequality, using $\ln(\lambda)/2\leq\tilde{r}\leq n^{2}\ln(\lambda)/(50(n-2k)^{2})$, we estimate

[TABLE]

Recalling $F(u+\delta)\geq\lceil\delta\rceil B(\lceil u+2\delta\rceil)$ and $\lceil\delta\rceil B(\lceil u\rceil)\geq F(u)-F(u+\delta)$, we compute

[TABLE]

Since $F(u)\geq 0.25$ and $\lambda>148$, we obtain

[TABLE]

Using $\Pr(X^{+}-X^{-}\geq 0.05\ln\lambda)\geq\lambda^{-0.8}/36$, we bound

[TABLE]

Finally $E(\Delta\mid k)\geq 0.02\cdot 0.05\ln\lambda=10^{-3}\ln\lambda$. ∎

We now extend the lemma to the whole region of $n/\ln\lambda\leq k<n/2$. If $k<2n/5$ the situation becomes easier because $4\leq c_{2}(k)<100$ and every $r$ in the smaller range $[c_{1}(k)\ln\lambda,\ln(\lambda)/2]$ provides at least an expected fitness improvement that is logarithmic in $\lambda$. Together with the previous lemma, we obtain the following statement for the drift in the whole region $n/\ln\lambda\leq k<n/2$.

Lemma 11**.**

Let $n$ be sufficiently large, $n/\ln\lambda\leq k<n/2$ and $\lambda\geq e^{5}$. If $c_{1}(k)\ln\lambda\leq r\leq c_{2}(k)\ln(\lambda)/100$ with $c_{1}(k),c_{2}(k)$ defined as in Lemma 9, then

[TABLE]

Proof.

If $r>\ln(\lambda)$, then $c_{2}(2n/5)=100$ implies $k>2n/5$ and the claim follows from Lemma 10. Hence let us assume $r\leq\ln(\lambda)$ in the remainder and compute $\tilde{\Delta}(\lambda/2,k,\tilde{r})$ with $\tilde{r}=r/2$.

We consider the probability $Q(k,i,\tilde{r})$ of creating from a parent with distance $k$ an offspring with fitness distance at most $k-i$ via standard bit mutation with mutation rate $p=\tilde{r}/n$. Let $i:=\max\{1,\lfloor c_{1}(k)\ln(\lambda)/2\rfloor\}\leq\max\{1,r/2\}=\tilde{r}$. If $\lfloor c_{1}(k)\ln(\lambda)/2\rfloor<1$, using $(1-p)^{n}\geq e^{-pn-p^{2}n}\geq(1-o(1))e^{-\tilde{r}}$ by Lemma 8(c), we compute

[TABLE]

Otherwise if $\lfloor c_{1}(k)\ln(\lambda)/2\rfloor>1$, using $\binom{k}{i}\geq(k/i)^{i}$, we obtain

[TABLE]

Hence, $\Pr(\tilde{\Delta}(\lambda/2,k,\tilde{r})\geq i\mid k)\geq 1-(1-1/\lambda)^{\lambda/2}\geq 1-e^{-1/2}>0.3$. We notice that $i\geq c_{1}(k)\ln(\lambda)/4$, since $i=1\geq c_{1}(k)\ln(\lambda)/4$ when $c_{1}(k)\ln\lambda<2$ and $i=\lfloor c_{1}(k)\ln(\lambda)/2\rfloor\geq c_{1}(k)\ln(\lambda)/4$ when $c_{1}(k)\ln\lambda\geq 2$. Consequently

[TABLE]

∎

If we only consider generations that use a rate within the right region, we can bound the expected runtime to reach $k\leq n/\ln\lambda$ by $O(n/\log\lambda)$ since the drift on the fitness is of order $\log\lambda$. The following theorem shows that the additional time spent with adjusting the rate towards the right region does not change this bound on the expected runtime.

Theorem 12**.**

The (1+$\lambda$) EA with self-adjusting mutation rate reaches a OneMax-value of $k\leq n/\ln\lambda$ within an expected number of $O(n/\log\lambda)$ iterations. This bound is valid regardless of the initial mutation rate.

Proof.

Let us denote by $k(x)$ the fitness distance of a search point $x\in\{0,1\}^{n}$.

We first argue that it takes an expected number of at most $O(\sqrt{n})$ iterations to reach a fitness distance of less than $n/2$. To this end, consider a parent $x$ with fitness distance $k(x)\geq n/2$. Let $2\leq r\leq n/4$ be the current rate and let $\tilde{r}\in\{r/2,2r\}$. Let $y$ be an offspring obtained from $x$ via standard bit mutation with mutation rate $\tilde{r}/n$. Let

[TABLE]

Then the fitness improvement is $k(x)-k(y)=X^{+}-X^{-}$ and both $X^{+}$ and $X^{-}$ follow binomial distributions. From elementary properties of the binomial distribution, see, e.g., [11], we have $\Pr(X^{+}\geq E(X^{+})+1)=\Omega(1)$ and $\Pr(X^{-}\leq E(X^{-}))=\Omega(1)$, and these are independent events. Hence with constant probability, we have

[TABLE]

Clearly, for the best offspring $z$ out of the $\lambda$ offspring generated in this iteration, we have $k(z)\leq k(y)$. Consequently, in each iteration starting with a parent $x$ with $k(x)\geq n/2$, with constant probability we gain at least one fitness level. By the additive drift theorem (see Theorem 3 and the subsequent discussion), from the initial random search point $x_{0}$ it takes an expected $O(\max\{0,k(x_{0})-n/2+1\})$ iterations to reach a parent with fitness distance less than $n/2$. Since $X:=k(x_{0})\sim\operatorname{Bin}(n,1/2)$, we have

[TABLE]

where we first exploited the symmetry of $X$ and then the well-known estimate $E(Y)^{2}\leq E(Y^{2})$ valid for all random variables $Y$, in particular, for $Y=|X-E(X)|$. By the law of total expectation, the expected time to reach a search point with $k$-value below $n/2$ is $O(\sqrt{n})$.

Without loss of generality we can now assume $k<n/2$ for the initial state. Our intuition is that once we begin to use a rate $r\in[(c_{1}(k)/2)\ln\lambda,c_{2}(k)\ln\lambda]$ at some distance level $k$, we will have a considerable drift on the OneMax-value and the strong drift on the rate keeps $r$ within or close to this interval. After we make progress and $k$ decreases to a new level, the corresponding $c_{1}$ and $c_{2}$ decrease, and the algorithm may take some time to readjust $r$ into new bounds.

We consider the stochastic process $X_{t}=\lfloor\log_{2}(r_{t})\rfloor$ and the current OneMax-value $Y_{t}$. According to Lemma 9 (a) and (b), there exists $\varepsilon=\Omega(1)$ such that

[TABLE]

Let $k_{0}>k_{1}>\dots>k_{N}$ be all the different OneMax-values taken by $Y_{t}$ until for the first time $Y_{t}\leq n/\ln\lambda$. By the additive drift theorem, it takes at most $O(\log n)$ iterations to have $r_{t}\in[(c_{1}(k_{0})/2)\ln\lambda,c_{2}(k_{0})\ln\lambda]$, regardless of how we set the initial rate. Since $c_{1}(Y_{t})$ is non-increasing, $r_{t}\geq c_{1}(Y_{t})\ln(\lambda)/2$ implies $r_{t}\geq c_{1}(Y_{t+1})\ln(\lambda)/2$. Thus, no readjustments are necessary to satisfy the lower bound $r_{t}\geq c_{1}(Y_{t})\ln(\lambda)/2$ once we have $r_{t}\geq c_{1}(k_{0})\ln(\lambda)/2$. We now regard the upper bound condition $r_{t}\leq c_{2}(Y_{t})\ln(\lambda)$. Let $(k_{s_{i}})$ be the subsequence consisting of all $k_{s_{i}}$ such that $\{t\mid r_{t}\leq c_{2}(Y_{t})\ln(\lambda),Y_{t}=k_{s_{i}}\}\neq\emptyset$. If $k_{i}$ is in the subsequence, once $r_{t}\leq c_{2}(Y_{t})\ln(\lambda)$ is achieved, the runtime until the first time $Y_{t}=k_{i+1}$ can be computed using the occupation lemma and the variable drift theorem. After that, it takes at most $O(\log(c_{2}(k_{i})/c_{2}(k_{i+1})))$ more iterations to readjust $r_{t}$ from $c_{2}(k_{i})\ln(\lambda)$ to $c_{2}(k_{i+1})\ln(\lambda)$ by the additive drift theorem. Similar arguments hold for $k_{i+1}$ if it is also in the subsequence. Otherwise let $k_{j}<k_{i+1}$ denote the next OneMax-values after $k_{i}$ in the subsequence. By definition, the fitness distance decreases from $k_{i+1}$ to $k_{j}$ before the rate being adjusted into the corresponding range. This means that improvements in distance are made during readjustment. Analogous to above, the expected (readjusting) time to decrease the rate from $c_{2}(k_{i+1})\ln(\lambda)$ to $c_{2}(k_{j})\ln(\lambda)$ is $O(\log(c_{2}(k_{i+1})/c_{2}(k_{j})))$. Therefore, the expected number of total readjusting time is at most

[TABLE]

When computing the expected number of non-adjusting iterations until $Y_{t}=k_{N}$, we choose a constant $b\in\mathbb{N}_{\geq 2}$ large enough such that $2e^{-2b\varepsilon/3}\leq 1/2-\delta/2$ holds for some positive constant $\delta>0$. When $r_{t}\in[(c_{1}(k)/2)\ln\lambda,c_{2}(k)\ln\lambda]$ and $Y_{t}=k$ for some $t$, then by Lemma 7, we obtain

[TABLE]

Hence, we obtain that $r_{t^{\prime}}\in[2^{-b-2}c_{1}(k)\ln\lambda,2^{b+1}c_{2}(k)\ln\lambda]$ happens with probability at least $\delta$. We see that there are at most $\lceil\log_{2}(100)\rceil$ steps between being in the range $[(c_{1}(k)/2)\ln\lambda,c_{2}(k)\ln\lambda]$ and being in the smaller range $[c_{1}(k)\ln\lambda,c_{2}(k)\ln(\lambda)/100]$ which is described in Lemma 11. If $r_{t}\in[2^{-b-2}c_{1}(k)\ln\lambda,2^{b+1}c_{2}(k)\ln\lambda]$, it takes at most a constant number of iterations $\alpha$ in expectation to reach $[c_{1}(k)\ln\lambda,c_{2}(k)\ln(\lambda)/100]$ because our mutation scheme employs a $50\%$ chance to perform a random step of the mutation rate. Based on Lemma 11, the fitness drift at distance $k$ of all rates within the narrow region is at least $10^{-3}\ln(\lambda)/\ln(en/k)$. This contributes to an average drift of at least

[TABLE]

for all random rates at distance $k$. Using the variable drift theorem (Theorem 3), we estimate the runtime as

[TABLE]

We notice that the expected runtime for adjusting $r_{t}$ is $O(\log n)=O(n/\log\lambda)$. Therefore, the total runtime is $O(n/\log\lambda)$ in expectation.

∎

6 Middle Region

In this section we estimate the expected number of generations until the number of one-bits has decreased from $k\leq n/\ln\lambda$ to $k\leq n/\lambda$. We first claim that the right region for $r$ is $1\leq r\leq\ln(\lambda)/2$. Hence, the (1+$\lambda$) EA is not very sensitive to the choice of $r$ here. Intuitively, this is due to the fact that a total fitness improvement of only $O(n/\log\lambda)$ suffices to cross the middle region, whereas an improvement of $\Omega(n)$ is needed for the far region.

We estimate the drift of the fitness in Lemma 13 and apply that result afterwards to estimate the number of generations to cross the region.

Lemma 13**.**

Let $n/\lambda\leq k\leq n/\ln\lambda$, $\lambda\geq 26$ and $2\leq r\leq\ln\lambda$. Then

[TABLE]

Proof.

We aim at computing $\tilde{\Delta}(\lambda/2,k,\tilde{r})$ with $\tilde{r}=r/2$. The probability that no zero-bit flips in a single offspring stemming from rate $\tilde{r}$ is $(1-\tilde{r}/n)^{n-k}\geq e^{-\tilde{r}}\geq 1/\sqrt{\lambda}$. We regard the number $Z$ of offspring of rate $\tilde{r}$ that have no flipped zeros. The expectation $E(Z)$ is at least $(1/\sqrt{\lambda})\cdot(\lambda/2)=\sqrt{\lambda}/2$. Applying Chernoff bounds (Theorem 5), we observe that $Z$ exceeds $\lambda_{0}:=\sqrt{\lambda}/4$ with probability at least $1-\exp(-\sqrt{\lambda}/16)>1/4$ since $\lambda\geq 26$. Assuming this to happen, we look at the first $\lambda_{0}$ offspring without flipped zeros. For $i\in\{1,\dots,\lambda_{0}\}$ let $X_{i}$ be the number of flipped ones in the $i$-th offspring. Then the $X_{i}$ are i.i.d. with $X_{i}\sim\operatorname{Bin}(k,\tilde{r}/n)$. Let $X^{*}\coloneqq\max\{X_{i}\}$. The expectation of $X^{*}$ is analyzed in Gießen and Witt [27, Lemma 4, part 2]; more precisely the following statement is shown for constant $\tilde{r}$:

[TABLE]

Since $\tilde{r}$ is not necessarily assumed constant here, we generalize the result by closely following the proof in [27]. Note that $X^{*}\geq 1$ if at least one of the offspring does not flip a zero-bit. Hence,

[TABLE]

where the second inequality used the assumption $\lambda_{0}k\tilde{r}/n\geq\alpha$, the third one $e^{x}\geq 1+x$ for $x\in\mathbb{R}$. and the fifth one the equivalent to the previous inequality $e^{-x}\leq 1/(1+x)$. Hence, Lemma 4, part 2 in [27] holds also for arbitrary $\tilde{r}$.

Setting $\alpha\coloneqq\lambda_{0}k\tilde{r}/n$, we now distinguish between two cases. If $\alpha\geq 1$, we obtain $E(X^{*})\geq\alpha/(1+\alpha)\geq 1/2$; otherwise we have $1+\alpha<2$, hence $E(X^{*})\geq\alpha/2=\lambda_{0}k\tilde{r}/(2n)$. Thus,

[TABLE]

Hence, using the law of total probability, we obtain the lower bound on the drift for the middle region. ∎

We now use our result on the drift to estimate the time spent in this region. We notice that $c_{2}(k)=4+o(1)$ when $k=o(n)$. This means we will often have $r_{t}\in[2,\ln(\lambda)]$ which provides the drift we need.

Theorem 14**.**

Let $\lambda\geq 26$ and $\lambda=n^{O(1)}$. Assume $k\leq n/\ln\lambda$ for the current OneMax-value of the self-adjusting (1+$\lambda$) EA. Then the expected number of generations until $k\leq n/\lambda$ is $O(n/\log\lambda)$.

Proof.

For $k\leq n/\ln\lambda$ the upper bound from Lemma 9 is $c_{2}(k)=4/(1-2/\ln\lambda)<11$. According to the lemma we have for $X_{t}\coloneqq\lceil\log_{2}r_{t}\rceil$ that $E(X_{t}-X_{t+1}\mid X_{t};X_{t}>\log_{2}(c_{2}(k)\ln\lambda))\geq 2\varepsilon=\Omega(1)$. The additive drift theorem yields that in $O(\log n)$ time we have $r_{t}\leq c_{2}(k)\ln\lambda$. We choose $b$ large enough such that $2e^{-2b\varepsilon/3}\leq 1-\delta$ holds for some positive constant $\delta>0$ and note that $b$ is constant. Applying Lemma 7 we obtain $\operatorname{Pr}(r_{t}\geq 2^{b}c_{2}(k)\ln\lambda)\leq 2e^{-2b\varepsilon/3}$ and $r_{t}\leq 2^{b}c_{2}(k)\ln\lambda$ happens with probability at least $\delta$. Once $r_{t}<2^{b}c_{2}(k)\ln\lambda<2^{b}(11\ln\lambda)$ it takes at most a constant number of iterations $\alpha$ in expectation to draw $r_{t}$ to $\ln(\lambda)$ or less. According to Lemma 13 this ensures a drift of at least $(1/4)\min\{1/2,\sqrt{\lambda}k/(8n)\}\geq(1/32)\min\{1,\sqrt{\lambda}k/n\}$, which implies an average drift of at least $c\min\{1,\sqrt{\lambda}k/n\}$ over all random rates at distance $k$, where $c>0$ is a constant. Considering $\min\{1,\sqrt{\lambda}k/n\}$, the minimum is taken on the first argument if $k>n/\sqrt{\lambda}$, and on the second if $k<n/\sqrt{\lambda}$.

We are interested in the expected time to reduce the OneMax-value to at most $n/\lambda$. To ease the application of drift analysis, we artificially modify the process and make it create the optimum when the state (OneMax-value) is strictly less than $n/\lambda$. Clearly, the first hitting time of state at most $n/\lambda$ does not change by this modification. Applying the variable drift theorem (Theorem 3) with $x_{\mathrm{min}}=n/\lambda$, $X_{0}=k\leq n/\ln\lambda$ and $h(x)=c\min\{1,\sqrt{\lambda}k/n\}$, the expected number of generations to reach state at most $n/\lambda$ is bounded from above by

[TABLE]

which is $O(n/\log\lambda)$. The overall expected number of generations spent is $O(\log n+n/\log\lambda)=O(n/\log\lambda)$ since $\lambda=n^{O(1)}$ by assumption. ∎

7 Near Region

In the near region, we have $k\leq n/\lambda$. Hence, the fitness is so low that we can expect only a constant number of offspring to flip at least one of the remaining one-bits. This assumes constant rate. However, higher rates are detrimental since they are more likely to destroy the zero-bits of the few individuals flipping one-bits. Hence, we expect the rate to drift towards constant values, as shown in the following lemma.

Lemma 15**.**

Let $k\leq n/\lambda$, $\lambda\geq 45$ and $4\leq r_{t}\leq\ln(\lambda)/4$. Then the probability that $r_{t+1}=r_{t}/2$ is at least $0.5099$.

Proof.

To prove the claim we exploit the fact that only few one-bits are flipped in both subpopulations. Using $r\coloneqq r_{t}$, we shall argue as follows. With sufficiently high (constant) probability, (i) the $2r$-subpopulation contains no individual strictly better than the parent, that is, with fitness less than $k$, and (ii) all $2r$-offspring with fitness at least $k$ are identical to the parent. Conditional on this, either the $r/2$-population contains individuals with fitness less than $k$ and the winning individual surely stems from this subpopulation, or the $r/2$-population contains no better offspring. In the latter case, we argue that there are many more individuals with fitness exactly $k$ in the $r/2$-population than in the $2r$-population, which gives a sufficiently high probability for taking the winning individual from this side (as it is chose uniformly at random from all offspring with fitness $k$).

Let $N_{r/2}$ and $N_{2r}$ be the number of offspring that did not flip any zero-bits using rate $r/2$ and $2r$, respectively. Then $E(N_{r/2})=(\lambda/2)(1-r/(2n))^{n-k}\geq(1-o(1))\lambda e^{-r/2}/2$, since $k\geq 1$ and

[TABLE]

where we used that $r\leq\ln\lambda/4$ and $\lambda=n^{O(1)}$, i. e. $\ln\lambda\leq c\ln n$ for some constant $c$. Using $k\leq n/\lambda$ we get $E(N_{2r})=(\lambda/2)(1-2r/n)^{n-k}\leq(\lambda/2)e^{-2r(1-1/\lambda)}$. In fact, we can discriminate $N_{r/2}$ and $N_{2r}$ by using Theorem 5 in the following way: we have

[TABLE]

for sufficiently large $n$, since $r\leq(\ln\lambda)/4$ and $\lambda\geq 45$. Similarly, we obtain

[TABLE]

Note that $e^{-2r(1-1/\lambda)}<e^{-r/2}/4$ holds for all $r\geq 4$ and $\lambda\geq 45$. Since the offspring are generated independently, the events $N_{r/2}>\lambda e^{-r/2}/4$ and $N_{2r}<\lambda e^{-2r(1-\frac{1}{\lambda})}$ happen together with probability at least $(1-0.175)\cdot(1-0.312)\geq 0.567\eqqcolon 1-p_{\text{err}_{1}}$. Conditioning on this and by using a union bound the probability $p_{\text{err}_{2}}$ that at least one of the $N_{2r}$ offspring that do not flip any zero-bits flips at least one one-bit can be upper bounded by

[TABLE]

using $k/n\leq 1/\lambda$ in the first and $r\geq 4$ and $\lambda\geq 45$ in the last inequality. By using a union bound, we find the probability $p_{\text{err}_{3}}$ that at least one $2r$-offspring flips exactly one one-bit and exactly one zero-bit to be at most

[TABLE]

for sufficiently large $n$, using $k/n\leq 1/\lambda$ for the first inequality. The third inequality is due to $\ln x-x\leq-(1-e^{-1})x<-(3/5)x$ for all $x>0$ (see Lemma 8.a) and the last inequality stems from $r\geq 4$. The second inequality follows from

[TABLE]

using again $r\leq\ln\lambda\leq c\ln n$ for some constant $c$. Let $M_{r}$ be the number of such offspring. Any other fitness-decreasing flip-combinations of zeroes and ones in the $2r$-subpopulation require an offspring to flip at least two one-bits. The probability that such an offspring is created is at most

[TABLE]

using $k\leq n/\lambda$ and $r\leq\ln\lambda/4$ and the fact that $(\ln^{2}x)/x$ is decreasing for $x\geq e^{2}$ and $\lambda\geq 45>e^{2}$.

The events $N_{r/2}>\lambda e^{-r/2}/4$, $N_{2r}<\lambda e^{-2r(1-\frac{1}{\lambda})}$, $M_{r}=0$, and the event that no fitness-decreasing offspring is created in the $2r$-subpopulation are sufficient to ensure that the best individual is either surely from the $r/2$-population or chosen uniformly at random from the $N_{r/2}+N_{2r}$ offspring. Conditioning on these events, we have that the probability that the best offspring is chosen from the $r/2$-population is at least

[TABLE]

Hence, using a union bound for the error probabilities, the unconditional probability is at least

[TABLE]

∎

We note that the restriction $r_{t}\geq 4$ in the lemma above is not strictly necessary. Also for smaller $r_{t}$, the probability that the winning individual is chosen from the $r_{t}/2$-population is by an additive constant larger than $1/2$. Showing this, however, would need additional proof arguments as for smaller $r_{t}$, the event that both subpopulations contain individuals with fitness $k-1$ becomes more likely. We avoid this additional technicality by only arguing for $r_{t}\geq 4$, which is enough since any constant $r_{t}$ is sufficient for the fitness drift we need (since we do not aim at making the leading constant precise).

In the following proof of the analysis of the near region, we use the above lemma (with quite some additional arguments) to argue that the $r$-value quickly reaches $4$ or less and from then on regularly returns to this region. This allows to argue that in the near region we have a speed-up of a factor of $\Theta(\lambda)$ compared to the (1+1) EA, since every offspring only has a probability of $O(1/\lambda)$ of making progress (see also [13, 27]).

Theorem 16**.**

Assume $k\leq n/\lambda$ for the current OneMax-value of the self-adjusting (1+$\lambda$) EA. Then the expected number of generations until the optimum is reached is $O(n\log(n)/\lambda+\log n)$.

Proof.

The aim is to estimate the OneMax-drift at the points in time (generations) $t$ where $r_{t}=O(1)$. To bound the expected number of generations until the mutation rate has entered this region, we basically consider the stochastic process $Z_{t}:=\max\{0,\lfloor\log_{2}(r_{t}/a)\rfloor\}$, where $a\coloneqq 4$, which is the lower bound on $r_{t}$ from Lemma 15. However, as we do not have proved a drift of $Z_{t}$ towards smaller values in the region $L\coloneqq(\ln\lambda)/4\leq r_{t}\leq 16\ln\lambda\eqqcolon U$ (where $16$ is an upper bound on $c_{2}(k)$ from Lemma 11), we use the potential function

[TABLE]

assuming that $L$ and $U$ have been rounded down and up to the closest power of $2$, respectively. We note that $X_{t}\geq\log_{2}(r_{t}/a)-1=\log_{2}(r_{t}/(2a))$ if $r_{t}\leq L$ and $X_{t}\geq(L/U)^{2}\log_{2}(r_{t}/(2a))$ in all three cases.

The potential function has a slope of $1$ for $a\leq r_{t}\leq L$. Lemma 15 gives us the drift $E(X_{t}-X_{t+1}\mid X_{t};a\leq r_{t}\leq L)\geq(0.5099-0.4901)/(2a)=\Omega(1)$. The function satisfies $X_{t}(Z_{t})-X_{t}(Z_{t}-1)\geq 4(X_{t}(Z_{t+1})-X_{t}(Z_{t}))$ if $L<r_{t}<U$, which corresponds to the region where the probability of decreasing $r_{t}$ by a factor of $1/2$ has only be bounded from below by $1/4$ due to the random steps. Still, $E(X_{t}-X_{t+1}\mid X_{t};L<r_{t}<U)\geq(1/4)4^{-\log_{2}(U/L)-1}-(3/4)4^{-\log_{2}(U/L)}=\Omega(1)$ in this region due to the concavity of the potential function. Finally, $E(X_{t}-X_{t+1}\mid X_{t};r_{t}\geq U)=4^{-\log_{2}(U/L)-1}(1/(2a))\Omega(1)=\Omega(1)$ by Lemma 11. Hence, altogether $E(X_{t}-X_{t+1}\mid X_{t};r_{t}\geq a)\geq\kappa$ for some constant $\kappa>0$. As $X_{0}=O(\log n)$, additive drift analysis yields an expected number of $O(\log n)$ generations until for the first time $X_{t}=0$ holds, corresponding to $r_{t}\leq a$. We denote this hitting time by $T$.

We now consider an arbitrary point of time $t\geq T$. The aim is to show a drift on the OneMax-value, depending on the current OneMax-value $Y_{t}$, which satisfies $Y_{t}\leq n/\lambda$ with probability $1$. To this end, we will use Lemma 7. We choose $b$ large enough such that $2e^{-2b\cdot\kappa/4}\leq 1-\delta$ holds for some positive constant $\delta>0$ and note that $b$ is constant. We consider two cases. If $X_{t}\leq b$, which happens with probability at least $\delta$ according to the lemma, then the bound $X_{t}\geq(L/U)^{2}(\log_{2}(r_{t}/(2a)))$ implies $r_{t}\leq 2^{(U/L)^{2}}2a2^{b}$. Hence, we have $r_{t}=O(1)$ in this case and obtain a probability of at least

[TABLE]

to improve the OneMax-value by $1$, using that $Y_{t}=O(n/\lambda)$ and pessimistically assuming a rate of $2r_{t}$ in all offspring. If $X_{t}>b$, we bound the improvement from below by [math]. Using the law of total probability, we obtain

[TABLE]

Now a multiplicative drift analysis with respect to the stochastic process on the $Y_{t}$, more precisely Theorem 4 using $\delta=\Theta(\lambda/n)$ and minimum state $1$, gives an expected number of $O((n/\lambda)\log Y_{0})=O(n\log(n)/\lambda)$ generations until the optimum is found. Together with the expected number $O(\log n)$ until the $r$-value becomes at most $a$, this proves the theorem. ∎

8 Putting Everything Together

In this section, we put together the analyses of the different regimes to prove our main result.

Proof of Theorem 1.

The lower bound actually holds for all unbiased parallel black-box algorithms, as shown in [2].

We add up the bounds on the expected number of generations spent in the three regimes, more precisely we add up the bounds from Theorem 12, Theorem 14 and Theorem 16, which gives us $O(n/\log\lambda+n\log(n)/\lambda+\log n)$ generations. Due to our assumption $\lambda=n^{O(1)}$ the bound is dominated by $O(n/\log\lambda+n\log(n)/\lambda)$ as suggested. ∎

Proof of Lemma 2.

We basically revisit the regions of different OneMax-values analyzed in this paper and bound the time spent in these regions under the assumption $r=\ln(\lambda)/2$. In the far region, Lemmas 10 and 11, applied with this value of $r$, imply a fitness drift of $\Omega(\log(\lambda)/\log(en/k))$ per generation, so the expected number of generations spent in the far region is $O(n/\log\lambda)$ as computed by variable drift analysis in the proof of Theorem 12.

The middle region is shortened at the lower end. For $k\geq n/\sqrt{\lambda}$, Lemma 13 gives a fitness drift of $\Omega(1)$, implying by additive drift analysis $O(n/\log\lambda)$ generations to reduce the fitness to at most $n/\sqrt{\lambda}$.

In the near region, which now starts at $n/\sqrt{\lambda}$, we have to argue slightly differently. Note that every offspring has a probability of at least $(1-r)^{n}\geq e^{-\ln(\lambda)/2+O(1)}=\Omega(\lambda^{-1/2})$ of not flipping a zero-bit. Hence, we expect $\Omega(\sqrt{\lambda})$ such offspring. We pessimistically assume that the other individuals do not yield a fitness improvement; conceptually, this reduces the population size to $\Omega(\sqrt{\lambda})$ offspring, all of which are guaranteed not to flip a zero-bit. Adapting the arguments from the proof of Theorem 16, the probability that at least of one of these individuals flips at least a one-bit is at least

[TABLE]

which is a lower bound on the fitness drift. Using the multiplicative drift analysis, the expected number of generations in the near region is $O(n\log(n)/\sqrt{\lambda})$. Putting the times for the regions together, we obtain the lemma. ∎

9 Experiments

Since our analysis is asymptotic in nature we performed some elementary experiments in order to see whether besides the asymptotic runtime improvement (showing an improvement for an unspecified large problem size $n$) we also see an improvement for realistic problem sizes. For this purpose we implemented the (1+$\lambda$) EA in C using the GNU Scientific Library (GSL) for the generation of pseudo-random numbers. In this section our performance measure is the runtime, represented by the number of generations until the optimum is found for the first time.

The first plot in Figure 1 displays the average runtime over 10000 runs of the self-adjusting (1+$\lambda$) EA on OneMax for $n=5000$ as given in Algorithm 1 over $\lambda=100,200,\ldots,1000$; the second plot shows the corresponding interquartile ranges to support that the results are statistically significant. We set the initial mutation rate to $2$, i. e., the minimum mutation rate the algorithm can attain. Moreover, the plot displays the average runtimes of the classic (1+$\lambda$) EA using static mutation probabilities of $1/n$ and of $(\ln\lambda)/(2n)$, the latter of which is asymptotically optimal for large $\lambda$ according to Lemma 2.

The average runtimes of both algorithms profit from higher offspring population sizes $\lambda$ leading to lower average runtimes as $\lambda$ increases. Interestingly, the two static settings of the classic (1+$\lambda$) EA outperform the self-adjusting (1+$\lambda$) EA for small values of $\lambda$ up to $\lambda=400$. For higher offspring population sizes the self-adjusting (1+$\lambda$) EA outperforms the classic ones, indicating that the theoretical performance gain of $\ln\ln\lambda$ can in fact be relevant in practice. Furthermore, we implemented the self-adjusting (1+$\lambda$) EA without the random steps, that is, when the rate is always adjusted according to how the best offspring are distributed over the two subpopulations. The experiments show that this variant of the self-adjusting (1+$\lambda$) EA performs generally slightly better on OneMax. Since the OneMax fitness landscape is structurally very simple, this result is not totally surprising. It seems very natural that the fitness of the best of $\lambda/2$ individuals, viewed as a function in the rate, is a unimodal function. In this case, the advantage of random steps to be able to leave local optima of this function is not needed. On the other hand, of course, this observation suggests to try to prove our performance bound rigorously also for the case without random rate adjustments. We currently do not see how to do this. Lastly, we implemented the (1+$\lambda$) EA using the fitness-depending mutation rate $p=\max\{\frac{\ln\lambda}{n\ln(en/d)},\frac{1}{n}\}$ as presented in [2]. The experiments suggest that this scheme outperforms all other variants considered.

Additionally we implemented another variant of the self-adjusting (1+$\lambda$) EA using three equally-sized subpopulations i. e. the additional one is using (that is, exploiting) the current mutation rate. We compared this variant with the self-adjusting (1+$\lambda$) EA, both with and without using random steps. The results are shown in Figure 2, again with respect to average runtimes and interquartile ranges. The experiments suggest that the variant using three subpopulations outperforms the self-adjusting (1+$\lambda$) EA slightly for small population sizes. For very high population sizes, using just two subpopulations seems to be a better choice.

To gain some understanding on how the parameters influence the runtime, we implemented the self-adjusting (1+$\lambda$) EA using different mutation rate update factors, that is, we consider the self-adjusting (1+$\lambda$) EA as given in Algorithm 1 where the mutation rate $r_{t}$ is increased or decreased by some factor $F$ (instead of the choice $F=2$ made in Algorithm 1). Note that we do not change the rule that we use the rates $r/2$ and $2r$ to create the subpopulations. Furthermore, after initialization, the algorithm starts with rate $F$ and the rate is capped below by $F$ and above by $1/(2F)$ during the run, accordingly.

The results are shown in Figure 3. The plot displays the average runtime over 10000 runs of the self-adjusting (1+$\lambda$) EA on OneMax for $n=5000$ over $\lambda=100,200,\ldots,1000$ using the update factors $F=2.0,1.5,1.01$. The plot suggests that lower values of $F$ yield a better performance. This result is not immediately obvious. Clearly, a large factor $F$ implies that the rate changes a lot from generation to generation (namely by a factor of $F$). These changes prevent the algorithm from using a very good rate for several iterations in a row. On the other hand, a small value for $F$ implies that it takes longer to adjust the rate to value that is far from the current one. We also performed experiments for $F>2$ and observed even worse runtimes than for $F=2.0$. The figure does not display the outcomes of these additional experiments since this would impair the readability of the diagram.

Finally, to illustrate the nontrivial development of the rate during a run of the algorithm we plotted the rate of three single runs of the self-adjusting (1+$\lambda$) EA using different factors $F$ over the fitness in Figure 4. Since the algorithm initialized with the rate $F$, the rate increases after initialization and decreases again with decreasing fitness-distance to the optimum. The plot suggests that for higher values of $F$ the rate is more unsteady due to the greater impact of the rate adjustments while smaller rate updates yield a more stable development of the rate. Interestingly, for all three values of $F$, the rates seem to correspond to the same rate after the initial increasing phase from $F$. Note that this illustration does not indicate the actual runtime. In fact, the specific runtimes are $19766$ for $F=1.01$, $19085$ for $F=1.05$ and $19857$ for $F=1.1$. A similar, more pronounced behaviour can be seen for $F=2.0$; we chose these particular values of $F$ for illustrative purposes since for $F=2.0$ the variance in the rate can be visually confusing for the reasons given above.

While we would draw from this experiment the conclusion that a smaller choice of $F$ is preferable in a practical application of our algorithm, the influence of the parameter on the runtime is not very large. So it might not be worth optimizing it and rather view Algorithm 1 as a parameter-less algorithm.

10 Conclusions

We proposed and analyzed a new simple self-adjusting mutation scheme for the (1+$\lambda$) EA. It consists of creating half the offspring with a slightly larger and the rest with a slightly smaller mutation rate. Based on the success of the subpopulations, the mutation rate is adjusted. This simple scheme overcomes difficulties of previous self-adjusting choices, e.g., the careful choice of the exploration-exploitation balance and the forgetting rate in the learning scheme of [20].

We proved rigorously that this self-adjusting (1+$\lambda$) EA optimizes the OneMax test function in an expected number of $O(n\lambda/\log\lambda+n\log n)$ fitness evaluations. This matches the runtime shown in [2] for a careful fitness-dependent choice of the mutation rate, which was also shown to be asymptotically optimal among all $\lambda$-parallel black-box optimization algorithms. Hence our runtime result indicates that the self-adjusting mechanism developed in this work is able to find very good mutation rates. To the best of our knowledge, this is the first time that a self-adjusting choice of the mutation rate speeds up a mutation-based algorithm on the OneMax test function by more than a constant factor.

The main technical challenge in this work is to analyze the quality of the best offspring. In contrast to most previous runtime analyses, where only the asymptotic order of the fitness gain was relevant, we needed a much higher degree of precision as we needed to make statements about in which sub-population the best offspring is, or, in case of multiple best offspring, how they are distributed over the two subpopulations. Note that the quality of the best offspring is not as strongly concentrated around its expectation as, e.g., the average quality.

As a side-result of our analyses, we have observed that using a fixed rate of $r=\ln(\lambda)/2$ gives the bound $O(n/\log\lambda+n\log(n)/\lambda^{1/2})$, which is also asymptotically optimal unless $\lambda$ is small. However, this setting is far off the usual constant choice of $r$. It is the first time that a significantly larger mutation rate was shown to be useful in a simple mutation-based algorithm for a simple fitness landscape. Previously, it was only observed that larger mutation rates can be helpful to leave local optima [23].

From this work, a number of open problems arise. A technical challenge is to prove that our algorithm also without the random rate adjustments performs well. This requires an even more precise analysis of the qualities of the offspring in the two sub-populations, for which we currently do not have the methods. From the view-point of understanding the mutation rate for population-based algorithms, two interesting questions are (i) to what extent our observation that larger mutation rates are beneficial for the (1+$\lambda$) EA on OneMax generalizes to other algorithms and problems, and (ii) for which other problems our self-adjusting choice of the mutation rate gives an improvement over the classic choice of $1/n$ or other static choices.

Acknowledgments

This work was supported by a public grant as part of the Investissement d’avenir project, reference ANR-11-LABX-0056-LMH, LabEx LMH, and by a grant by the Danish Council for Independent Research (DFF-FNU 4002–00542).