Controlling the Stopping Criteria: BI-POP CMA-ES

A variant of the Covariance Matrix Adaptation Evolution Strategy (CMA-ES) [Hansen2001] implies to control very specifically the termination criteria in the generational loop. This can be achieved by writing the algorithm partially invoking manually the generate() and update() inside a loop with specific stop criteria. In fact, the BI-POP CMA-ES [Hansen2009] has 9 different stop criteria, which are used to control the independent restarts, with different population sizes, of a standard CMA-ES.

As usual, the first thing to do is to create the types and as usual, we’ll need a minimizing fitness and an individual that is a list.

N = 30

The main function includes the setting of some parameters, namely the number of increasing population restarts and the initial sigma value. Then, the instanciation of the Toolbox is done in the main function because it will change with the restarts. Next are initialized the HallOfFame, The statistics and the list of Logbook objects, one for each restart.

creator.create("Individual", list, fitness=creator.FitnessMin)

def main(verbose=True):
    NRESTARTS = 10  # Initialization + 9 I-POP restarts
    SIGMA0 = 2.0    # 1/5th of the domain [-5 5]

    toolbox = base.Toolbox()
    toolbox.register("evaluate", benchmarks.rastrigin)

    halloffame = tools.HallOfFame(1)
    stats = tools.Statistics(lambda ind: ind.fitness.values)
    stats.register("avg", numpy.mean)
    stats.register("std", numpy.std)
    stats.register("min", numpy.min)
    stats.register("max", numpy.max)
    
    logbooks = list()

Then the first loop controlling the restart is set up. It encapsulates the generational loop with its many stop criteria. The content of this last loop is simply the generate-update loop as presented in the deap.algorithms.eaGenerateUpdate() function.

    while i < (NRESTARTS + nsmallpopruns):
        strategy = cma.Strategy(centroid=numpy.random.uniform(-4, 4, N), sigma=sigma, lambda_=lambda_)
        toolbox.register("generate", strategy.generate, creator.Individual)
        toolbox.register("update", strategy.update)
        
        logbooks.append(tools.Logbook())
        logbooks[-1].header = "gen", "evals", "restart", "regime", "std", "min", "avg", "max"
        
        conditions = {"MaxIter" : False, "TolHistFun" : False, "EqualFunVals" : False,
                      "TolX" : False, "TolUpSigma" : False, "Stagnation" : False,
                      "ConditionCov" : False, "NoEffectAxis" : False, "NoEffectCoor" : False}
        while not any(conditions.values()):
            # Generate a new population
            population = toolbox.generate()
            
            # Evaluate the individuals
            fitnesses = toolbox.map(toolbox.evaluate, population)
            for ind, fit in zip(population, fitnesses):
                ind.fitness.values = fit
            
            halloffame.update(population)
            record = stats.compile(population)
            logbooks[-1].record(gen=t, evals=lambda_, restart=i, regime=regime, **record)
            if verbose:
                print(logbooks[-1].stream)

            # Update the strategy with the evaluated individuals
            toolbox.update(population)
            if t >= MAXITER:
                # The maximum number of iteration per CMA-ES ran
                conditions["MaxIter"] = True
            
            mins.append(record["min"])
            if (len(mins) == mins.maxlen) and max(mins) - min(mins) < TOLHISTFUN:
                # The range of the best values is smaller than the threshold
                conditions["TolHistFun"] = True

            if t > N and sum(equalfunvalues[-N:]) / float(N) > EQUALFUNVALS:
                # In 1/3rd of the last N iterations the best and k'th best solutions are equal
                conditions["EqualFunVals"] = True

            if all(strategy.pc < TOLX) and all(numpy.sqrt(numpy.diag(strategy.C)) < TOLX):
                # All components of pc and sqrt(diag(C)) are smaller than the threshold
                conditions["TolX"] = True
            
            # Need to transfor strategy.diagD[-1]**2 from pyp/numpy.float64 to python
            # float to avoid OverflowError
            if strategy.sigma / sigma > float(strategy.diagD[-1]**2) * TOLUPSIGMA:
                # The sigma ratio is bigger than a threshold
                conditions["TolUpSigma"] = True
            
            if len(bestvalues) > STAGNATION_ITER and len(medianvalues) > STAGNATION_ITER and \
               numpy.median(bestvalues[-20:]) >= numpy.median(bestvalues[-STAGNATION_ITER:-STAGNATION_ITER + 20]) and \
               numpy.median(medianvalues[-20:]) >= numpy.median(medianvalues[-STAGNATION_ITER:-STAGNATION_ITER + 20]):
                # Stagnation occured
                conditions["Stagnation"] = True

            if strategy.cond > 10**14:
                # The condition number is bigger than a threshold
                conditions["ConditionCov"] = True

            if all(strategy.centroid == strategy.centroid + 0.1 * strategy.sigma * strategy.diagD[-NOEFFECTAXIS_INDEX] * strategy.B[-NOEFFECTAXIS_INDEX]):
                # The coordinate axis std is too low
                conditions["NoEffectAxis"] = True

            if any(strategy.centroid == strategy.centroid + 0.2 * strategy.sigma * numpy.diag(strategy.C)):
        stop_causes = [k for k, v in conditions.items() if v]
        print("Stopped because of condition%s %s" % ((":" if len(stop_causes) == 1 else "s:"), ",".join(stop_causes)))
        i += 1

Some variables have been omited for clarity, refer to the complete example for more details examples/es/cma_bipop.

[Hansen2001]Hansen and Ostermeier, 2001. Completely Derandomized Self-Adaptation in Evolution Strategies. Evolutionary Computation
[Hansen2009]Hansen, 2009. Benchmarking a BI-Population CMA-ES on the BBOB-2009 Function Testbed. GECCO‘09