Parameter Control and Adaptation.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "# Parameter Control and Adaptation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "outputs": [],
   "source": [
    "import random\n",
    "import pprint\n",
    "import math\n",
    "import copy\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# For presenting as slides\n",
    "#plt.rcParams['figure.figsize'] = [12, 8]\n",
    "#plt.rcParams.update({'font.size': 22})\n",
    "#plt.rcParams['lines.linewidth'] = 3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Evolutionary algorithms have many different parameters and implementation choices. There are qualitative parameters in the design of an algorithm (e.g., what representation, what recombination operators, what selection operators, etc.) and quantitative parameters (e.g., mutation rate, crossover rate, selection bias, population size, etc.) So far, we used intuition and common default values for many of these parameters. However, parameters can have a large influence on the performance of an algorithm on a particular problem. This chapter looks at how to optimise the selection of parameter values."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Example Problem: Pairwise Testing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "As an example search problem we will consider pairwise testing, a form of combinatorial interaction testing where we want to generate covering arrays of strength 2. The problem is as follows: We are testing a system with _x_ parameters, each parameter can have a different number of parameter values. In order to properly test the system, we would need to check all possible combinations of parameter values, but the combinatorial explosion usually makes this practically impossible. Since most failures are triggered by combinations of only few parameter values, the idea of pairwise testing is to try to cover all combinations of values for any pair of 2 parameters, which greatly reduces the number of tests."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "We'll use a very basic genetic algorithm to solve the problem. We'll use our wrapper class because we will be caching some values later on:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "class L(list):\n",
    "    \"\"\"\n",
    "    A subclass of list that can accept additional attributes.\n",
    "    Should be able to be used just like a regular list.\n",
    "    \"\"\"\n",
    "    def __new__(self, *args, **kwargs):\n",
    "        return super(L, self).__new__(self, args, kwargs)\n",
    "\n",
    "    def __init__(self, *args, **kwargs):\n",
    "        if len(args) == 1 and hasattr(args[0], '__iter__'):\n",
    "            list.__init__(self, args[0])\n",
    "        else:\n",
    "            list.__init__(self, args)\n",
    "        self.__dict__.update(kwargs)\n",
    "\n",
    "    def __call__(self, **kwargs):\n",
    "        self.__dict__.update(kwargs)\n",
    "        return self"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Let's assume we have a system with a number of parameters, each with a number of possible values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "parameters = { 0 : [\"a\", \"b\", \"c\", \"d\", \"x\"],\n",
    "               1 : [\"d\", \"e\", \"f\", \"g\", \"x\"],\n",
    "               2 : [\"g\", \"h\", \"i\", \"j\", \"x\"],\n",
    "               3 : [\"g\", \"h\", \"i\", \"j\", \"x\"],\n",
    "               4 : [\"g\", \"h\", \"i\", \"j\", \"x\"],\n",
    "               5 : [\"j\", \"k\", \"l\", \"m\", \"x\"]}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "A test case is a list of values, one for each parameter of our system:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def get_random_test():\n",
    "    num_params = len(parameters)\n",
    "    test = []\n",
    "    for p in range(num_params):\n",
    "        test.append(random.choice(parameters[p]))\n",
    "    return test"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "The aim of our optimisation problem is to produce a set of tests that covers as many as possible combinations of parameter-value pairs, at the same time we don't want too many test cases. In principle we could make this a multi-objective problem and optimise for size and coverage, but to simplify our examples we will just assume a fixed number of tests, and try to cover as many as possible pairs with that exact number of tests:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "num_tests = 25"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "An individual of our search is a list of `num_tests` of such tests:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def get_random_solution():\n",
    "    solution = L([])\n",
    "    for test in range(num_tests):\n",
    "        solution.append(get_random_test())\n",
    "    return solution"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "As a simple (though maybe somewhat computationally inefficient) fitness function we count how many distinct pairs of pairwise parameter-value combinations we are covering by putting them in a set and then counting the size. The larger the set, the more pairs we have covered:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def get_fitness(individual):\n",
    "    pairs = set()\n",
    "\n",
    "    for num1 in range(len(parameters) - 1):\n",
    "        for num2 in range(num1 + 1, len(parameters)):\n",
    "            for row in individual:\n",
    "                pairs.add((num1, num2, row[num1], row[num2]))\n",
    "\n",
    "    return len(pairs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "In previous chapters, the parameters were dispersed throughout the notebooks with some global variables. Since we are focusing on parameters now, we will store them explicitly in a dictionary. Let's set some default values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "configuration = {\n",
    "    \"P_xover\": 0.7,\n",
    "    \"P_mutation\": 0.07,\n",
    "    \"population_size\": 100,\n",
    "    \"tournament_size\": 2\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "We are only optimising quantitative parameters, so the operators will be the ones we've used previously, adapted only to make use of our configuration dictionary:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def tournament_selection(population):\n",
    "    # Make sure the sample isn't larger than the population\n",
    "    candidates = random.sample(population, min(len(population), configuration[\"tournament_size\"]))\n",
    "    winner = max(candidates, key=lambda x: get_fitness(x))    \n",
    "                \n",
    "    return winner"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "We'll use single-point crossover for tests, and just need to make sure to preserve our wrapper classes:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def crossover(parent1, parent2):\n",
    "    pos = random.randint(1, len(parent1))\n",
    "\n",
    "    offspring1 = L(copy.deepcopy(parent1[:pos] + parent2[pos:]))\n",
    "    offspring2 = L(copy.deepcopy(parent2[:pos] + parent1[pos:]))\n",
    "\n",
    "    return offspring1, offspring2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "In the past we dynamically set the mutation rate based on the size of individuals. To make it easier to evaluate the effects of the mutation rate parameter, we will now set it explicitly. Each parameter in our list of tests is replaced with a different parameter value with that probability. We'll pass in the probability as a parameter rather and read the value of the configuration dictionary outside the operator."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def mutate(individual, P_mutate):\n",
    "    mutated = L(copy.deepcopy(individual))\n",
    "    for num_row in range(len(mutated)):\n",
    "        for num_col in range(len(parameters)):\n",
    "            if random.random() < P_mutate:\n",
    "                choice = parameters[num_col][:]\n",
    "                choice.remove(mutated[num_row][num_col])\n",
    "                mutated[num_row][num_col] = random.choice(choice)\n",
    "\n",
    "    return mutated"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "The reproduction-related probabilities are considered during a step of our algorithm:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def evolution_step(population):\n",
    "    new_population = []\n",
    "    \n",
    "    while len(new_population) < len(population):\n",
    "        parent1 = selection(population)\n",
    "        parent2 = selection(population)\n",
    "\n",
    "        if random.random() < configuration[\"P_xover\"]:\n",
    "            offspring1, offspring2 = crossover(parent1, parent2)\n",
    "        else:\n",
    "            offspring1, offspring2 = parent1, parent2\n",
    "\n",
    "        offspring1 = mutate(offspring1, configuration[\"P_mutation\"])\n",
    "        offspring2 = mutate(offspring2, configuration[\"P_mutation\"])\n",
    "\n",
    "        new_population.append(offspring1)\n",
    "        new_population.append(offspring2)\n",
    "\n",
    "    population.clear()\n",
    "    population.extend(new_population)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Finally, we'll just set some parameters and global variables for our experiments, such as the list of fitness values to observe what happened, and the number of fitness evaluations as our stopping criterion:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "max_evaluations = 10000\n",
    "selection = tournament_selection\n",
    "fitness_values = []"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "The algorithm itself should by now be well known:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def ga():\n",
    "    population = [get_random_solution() for _ in range(configuration[\"population_size\"])]\n",
    "    best_individual = max(population, key=lambda k: get_fitness(k))\n",
    "    best_fitness = get_fitness(best_individual)\n",
    "    print(f\"Initial population, best fitness: {best_fitness}\")\n",
    "\n",
    "    iteration = 0\n",
    "    while iteration < (max_evaluations / configuration[\"population_size\"]):\n",
    "        fitness_values.append(best_fitness)\n",
    "        iteration += 1\n",
    "        evolution_step(population)\n",
    "        current_best = max(population, key=lambda k: get_fitness(k))\n",
    "        current_fitness = get_fitness(current_best)\n",
    "        if current_fitness > best_fitness:\n",
    "            print(f\"Iteration {iteration}, best fitness: {best_fitness}\")\n",
    "            best_individual = copy.deepcopy(current_best)\n",
    "            best_fitness = current_fitness\n",
    "\n",
    "    fitness_values.append(best_fitness)\n",
    "    return best_individual"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "We can now run our algorithm and look at the resulting test suite:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true,
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "fitness_values = []\n",
    "solution = ga()\n",
    "\n",
    "# Pretty print output\n",
    "pp = pprint.PrettyPrinter(depth=6)\n",
    "solution.sort()\n",
    "pp.pprint(solution)\n",
    "print(len(solution))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "plt.plot(fitness_values, label=\"Fitness\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Parameter Tuning"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Tuning refers to a (systematic) process of optimising the values of the parameters of our algorithm. There are many considerations when doing so: Which parameters, and which levels? We can't look at all possible values for our parameters, so we will have to consider some relevant values. For example, let's consider the following parameters and their possible values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "search_parameters = {\n",
    "    \"population_size\": [5, 37, 52, 67, 99],\n",
    "    \"tournament_size\": [1, 2, 3, 4, 5],\n",
    "    \"P_mutation\": [0.01, 0.34, 0.5, 0.65, 0.99],\n",
    "    \"P_xover\": [0.01, 0.34, 0.5, 0.65, 0.99]\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "To find the best configuration based on these options, one way is to produce a factorial design, i.e., a combination of all possible values (very similar to the example problem we are solving with search in this chapter!) Let's define a function that produces such a factorial design:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def get_factorial_design(search_parameters):\n",
    "    configurations = [{}]\n",
    "    \n",
    "    for parameter in search_parameters.keys():\n",
    "        new_configurations = []\n",
    "        for configuration in configurations:\n",
    "            for option in range(len(search_parameters[parameter])):\n",
    "                cp = configuration.copy()\n",
    "                cp[parameter] = search_parameters[parameter][option]\n",
    "                new_configurations.append(cp)\n",
    "            \n",
    "        configurations = new_configurations\n",
    "    \n",
    "    return configurations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "For our choice of parameters and values, this is the resulting factorial design:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "# This produces a loooong list...\n",
    "# get_factorial_design(search_parameters)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "len(get_factorial_design(search_parameters))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "In order to find the best parameter setting, we would have to run the search on some example problems for each of these configurations. That would take very long, so one option is to produce a fractional factorial design where we omit or restrict certain combinations. Let's simply reduce the value choices to 2 per parameter:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "search_parameters = {\n",
    "    \"population_size\": [20, 100],\n",
    "    \"tournament_size\": [1, 3],\n",
    "    \"P_mutation\": [0.01, 0.5],\n",
    "    \"P_xover\": [0.0, 0.7]\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "len(get_factorial_design(search_parameters))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "In order to evaluate how well each of the configurations performs, we need to define what _performance_ actually means. A simple measurement of performance is given by the fitness value of the resulting solution. However, since the algorithm is randomised, we need to apply it repeatedly. Furthermore, we need to choose which problem instances to evaluate the algorithm on. To keep things simple, we will only consider the specific problem listed earlier, and we'll measure the performance as the average of the best fitness values for a number of repetitions:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "from IPython.utils import io\n",
    "def get_performance(config):\n",
    "    global max_evaluations\n",
    "    results = []\n",
    "    repetitions = 4 # Should be more, but to accelerate the notebook...\n",
    "    old_evaluations = max_evaluations\n",
    "    max_evaluations = 3000 # Should be more, but to accelerate the notebook...\n",
    "    for i in range(repetitions):\n",
    "        fitness_values = []\n",
    "        configuration.update(config)\n",
    "        with io.capture_output() as captured: \n",
    "            result = ga()\n",
    "            fitness = get_fitness(result)\n",
    "        results.append(fitness)\n",
    "    max_evaluations = old_evaluations\n",
    "    return sum(results)/len(results)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Now we just have to consider each configuration of our factorial design and measure the performance. The best configuration is the one with the best performance:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "best_config = None\n",
    "best_performance = 0\n",
    "#max_evaluations = 10000 # Limit number of fitness evaluations so we don't have to wait too long\n",
    "for config in get_factorial_design(search_parameters):\n",
    "    print(f\"Current configuration: {config}\")\n",
    "    performance = get_performance(config)\n",
    "    print(f\"Performance: {performance}\")\n",
    "    if performance > best_performance:\n",
    "        best_config = config\n",
    "        best_performance = performance\n",
    "\n",
    "print(f\"Chosen best configuration: {best_config}\")\n",
    "configuration.update(best_config)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Finding vectors of parameter values is a complex optimisation task. Our performance function defines a _utility landscape_ -- an abstract landscape where the locations are the parameter vectors of an EA and the height reflects the utility (performance). Thus, technically we can treat our tuning problem as a search problem. To do so, we'll define range limits for the parameters we are optimising:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "search_parameter_ranges = {\n",
    "    \"population_size\": [1, 100],\n",
    "    \"tournament_size\": [1, 5],\n",
    "    \"P_mutation\": [0.0, 0.2],\n",
    "    \"P_xover\": [0.0, 1.0]\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "In order to apply meta-heuristic search to the problem of optimising the parameters of a meta-heuristic search algorithm, we need to define a representation, which is a parameter vector. We'll define a helper function that gives us a random value for a chosen parameter, using the range limits defined above:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def get_random_parameter(parameter_name):\n",
    "    value = random.uniform(*search_parameter_ranges[parameter_name])\n",
    "    if type(1) == type(search_parameter_ranges[parameter_name][0]):\n",
    "        value = int(value)\n",
    "    return value"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "An individual of our parameter-search is a dictionary of parameter values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def get_random_configuration():\n",
    "    config = {}\n",
    "    for param in search_parameter_ranges.keys():\n",
    "        config[param] = get_random_parameter(param)\n",
    "    return config"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "get_random_configuration()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "A very basic mutation of parameter dictionaries is to replace individual parameters with different values."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def mutate_configuration(config):\n",
    "    copy = config.copy()\n",
    "    P_mutate = 1/len(copy)\n",
    "    while copy == config:\n",
    "        for key in copy.keys():\n",
    "            if random.random() < P_mutate:\n",
    "                copy[key] = get_random_parameter(key)\n",
    "\n",
    "    return copy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Technically, we can use any search algorithm for this optimisation problem. However, measuring fitness requires running our `get_performance` function, which is expensive. Thus, some algorithms are better suited than others. One popular algorithm in practice is _differential evolution_, which we haven't covered yet, so to keep things simple we'll just use a (1+1)EA:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def oneplusone():\n",
    "    current = best_config #get_random_configuration()\n",
    "    fitness = get_performance(current)\n",
    "    iteration = 0\n",
    "    print(f\"Iteration {iteration}: Fitness {fitness} - {current}\")\n",
    "\n",
    "    while iteration < max_steps:\n",
    "        iteration += 1\n",
    "        candidate = mutate_configuration(current)\n",
    "        candidate_fitness = get_performance(candidate)\n",
    "        print(f\"Iteration {iteration}: Fitness {candidate_fitness} - {candidate}\")\n",
    "        if candidate_fitness >= fitness:\n",
    "            current = candidate\n",
    "            fitness = candidate_fitness\n",
    "        fitness_values.append(fitness)\n",
    "\n",
    "    return current"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true,
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "max_steps = 20\n",
    "tuned_opo = oneplusone()\n",
    "tuned_opo"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "To reduce the enormous costs of fitness evaluations in this search, the _response surface methodology_ systematically samples the utility space, creates a regression model that predicts the utility for parameter values, and then uses this regression model as a surrogate for querying the performance."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "As an alternative to a full factorial design, the response surface methodology is often performed using a _central composite design_. For the central composite design we obtain five levels for each parameter: {0, 1, −1, α, −α}.\n",
    "\n",
    "The design consists of three matrices:\n",
    "- A factorial design with two levels for each factor (1, -1)\n",
    "- A set of center points (0)\n",
    "- A set of axial points, configurations identical to the centre points except for one factor that is varied with values below and above the median (α, −α).\n",
    "\n",
    "(The choice of 1, -1, α, -α is actually more complicated, but out of scope for this notebook)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "For example, assume we have derived the following levels (-α, -1, 0, 1, α) for our factors:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "search_parameters = {\n",
    "    \"population_size\": [5, 37, 52, 67, 99],\n",
    "    \"tournament_size\": [1, 2, 3, 4, 5],\n",
    "    \"P_mutation\": [0.005, 0.01, 0.05, 0.1, 0.2],\n",
    "    \"P_xover\": [0.01, 0.34, 0.5, 0.65, 0.99]\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "We will define helper functions that produce the three component matrices:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def get_factorial_for(search_parameters, options = [1, 3]):\n",
    "    configurations = [{}]\n",
    "    \n",
    "    for parameter in search_parameters.keys():\n",
    "        new_configurations = []\n",
    "        for configuration in configurations:\n",
    "            for option in options:\n",
    "                cp = configuration.copy()\n",
    "                cp[parameter] = search_parameters[parameter][option]\n",
    "                new_configurations.append(cp)\n",
    "            \n",
    "        configurations = new_configurations\n",
    "    \n",
    "    return configurations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "get_factorial_for(search_parameters)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def get_center_point(parameters, center = 2):\n",
    "    configuration = {}\n",
    "    for parameter in parameters.keys():\n",
    "        configuration[parameter] = parameters[parameter][center]\n",
    "    return configuration"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "get_center_point(search_parameters)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def get_axial_points(parameters, axial_points = [0, 4]):\n",
    "    configurations = []\n",
    "    \n",
    "    for parameter in parameters.keys():\n",
    "        for axial_point in axial_points:\n",
    "            configuration = get_center_point(parameters)\n",
    "            configuration[parameter] = parameters[parameter][axial_point]\n",
    "            configurations.append(configuration)\n",
    "    return configurations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "get_axial_points(search_parameters)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "The full central composite design now consists of all three of these:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "configurations = get_factorial_for(search_parameters) + \\\n",
    "    [ get_center_point(search_parameters) ] + \\\n",
    "    get_axial_points(search_parameters)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "len(configurations)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Now we need some patience while we determine the performance for each of these configurations:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true,
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "responses = []\n",
    "max_evaluations = 10000\n",
    "for config in configurations:\n",
    "    r = get_performance(config)\n",
    "    config[\"performance\"] = r\n",
    "    print(config)\n",
    "    responses.append(r)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "The next step is to create a regression model that predicts the performance of any given combination of parameter values. We'll cut this short and simply throw our data at a multilayer perceptron regression model, and doing this properly and checking the fit of the model are left as an exercise."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "from sklearn.neural_network import MLPRegressor"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "model = MLPRegressor(solver=\"lbfgs\", max_iter=1000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "df = pd.DataFrame(configurations)\n",
    "y = df['performance'].values\n",
    "x = df[['population_size', 'tournament_size', 'P_mutation', 'P_xover']].values\n",
    "# This will need some adjusting...\n",
    "model.fit(x, y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "We can now take any configuration and query the model. For example, let's query the configuration we determined best with our simplified factorial design initially:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "model.predict([list(best_config.values())])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "We'll use a (1+1)EA again to do the search, but instead of querying the fitness function, we query the model:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "def surrogate_oneplusone():\n",
    "    current = best_config #get_random_configuration()\n",
    "    fitness = model.predict([list(current.values())])[0]\n",
    "    iteration = 0\n",
    "    print(f\"Iteration {iteration}: Fitness {fitness} - {current}\")\n",
    "\n",
    "    while iteration < max_steps:\n",
    "        iteration += 1\n",
    "        candidate = mutate_configuration(current)\n",
    "        candidate_fitness = model.predict([list(candidate.values())])[0]\n",
    "        if candidate_fitness >= fitness:\n",
    "            if candidate_fitness > fitness:\n",
    "                print(f\"Iteration {iteration}: Fitness {candidate_fitness} - {candidate}\")\n",
    "                \n",
    "            current = candidate\n",
    "            fitness = candidate_fitness\n",
    "        fitness_values.append(fitness)\n",
    "\n",
    "    return current"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Since querying the model is very cheap compared to measuring the actual performance, we can run many more iterations of the search algorithm."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "max_steps = 1000\n",
    "rsm_config = surrogate_oneplusone()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "rsm_config"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Is the optimised configuration really better? Let's find out by comparing against our initial default values."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "configuration['P_xover'] = 0.7\n",
    "configuration['P_mutation'] = 0.07\n",
    "configuration['population_size'] = 100\n",
    "configuration['tournament_size'] = 2\n",
    "\n",
    "fitness_values = []\n",
    "ga()\n",
    "# Stretch this by the population size to make plots comparable\n",
    "default_fitness = [item for item in fitness_values for i in range(configuration['population_size'])]\n",
    "\n",
    "\n",
    "configuration.update(tuned_opo)\n",
    "fitness_values = []\n",
    "ga()\n",
    "# Stretch this by the population size to make plots comparable\n",
    "opo_tuned_fitness = [item for item in fitness_values for i in range(configuration['population_size'])]\n",
    "\n",
    "\n",
    "configuration.update(rsm_config)\n",
    "fitness_values = []\n",
    "ga()\n",
    "# Stretch this by the population size to make plots comparable\n",
    "rsm_tuned_fitness = [item for item in fitness_values for i in range(configuration['population_size'])]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "plt.plot(opo_tuned_fitness, label=\"Tuned (1+1)\")\n",
    "plt.plot(rsm_tuned_fitness, label=\"Tuned (RSM)\")\n",
    "plt.plot(default_fitness, label=\"Default\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Parameter Control"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "There's actually something fundamentally missing in the tuning experiments: We did all our tuning on the exact same problem. That's obviously not what we would do in practice: In order to find good parameters we performed search on that problem; once we have found a good set of parameters we know that it is a good set of parameters because we have found a good solution already, so we don't need to run the search again on that problem. What is missing in our experiments is that we would in practice of course do the tuning on several problems, and then hope that these parameters generalise to other, new problems."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "A further limitation of parameter tuning is that it tries to select a good set of parameter values, but these parameters remain constant throughout the search. However, different phases of the search may benefit from different parameter values. For example, during initial phases of the search we may want to put more focus on exploration, thus using higher mutation rates than in later phases (like we saw with the MIO algorithm). Thus, a basic approach to control parameters is to make the parameter values dependent on the phase of the search:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def calculate_probability(size, iteration, max_iterations):\n",
    "    return math.pow(4 + (size - 2)/(max_iterations - 1) * iteration, -1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "mutation_values = [calculate_probability(num_tests*len(parameters), x, 100) for x in range(100)]\n",
    "plt.plot(mutation_values)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "We'll now create a version of the genetic algorithm which doesn't use the mutation rate from our configuration dictionary, but derives it using our probability function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def controlled_evolution_step(population, P_mutation):\n",
    "    new_population = []\n",
    "    \n",
    "    while len(new_population) < len(population):\n",
    "        parent1 = selection(population)\n",
    "        parent2 = selection(population)\n",
    "\n",
    "        if random.random() < configuration[\"P_xover\"]:\n",
    "            offspring1, offspring2 = crossover(parent1, parent2)\n",
    "        else:\n",
    "            offspring1, offspring2 = parent1, parent2\n",
    "\n",
    "        offspring1 = mutate(offspring1, P_mutation)\n",
    "        offspring2 = mutate(offspring2, P_mutation)\n",
    "\n",
    "        new_population.append(offspring1)\n",
    "        new_population.append(offspring2)\n",
    "\n",
    "    population.clear()\n",
    "    population.extend(new_population)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "mutation_values = []\n",
    "def controlled_ga():\n",
    "    population = [get_random_solution() for _ in range(configuration[\"population_size\"])]\n",
    "    best_individual = max(population, key=lambda k: get_fitness(k))\n",
    "    best_fitness = get_fitness(best_individual)\n",
    "    print(f\"Initial population, best fitness: {best_fitness}\")\n",
    "\n",
    "    iteration = 0\n",
    "    max_iterations = (max_evaluations / configuration[\"population_size\"])\n",
    "    while iteration < max_iterations:\n",
    "        fitness_values.append(best_fitness)\n",
    "        P_mutation = calculate_probability(num_tests * len(parameters), iteration, max_iterations)\n",
    "        mutation_values.append(P_mutation)\n",
    "        iteration += 1\n",
    "        controlled_evolution_step(population, P_mutation)\n",
    "        current_best = max(population, key=lambda k: get_fitness(k))\n",
    "        current_fitness = get_fitness(current_best)\n",
    "        if current_fitness > best_fitness:\n",
    "            print(f\"Iteration {iteration}, best fitness: {best_fitness}\")\n",
    "            best_individual = copy.deepcopy(current_best)\n",
    "            best_fitness = current_fitness\n",
    "\n",
    "    fitness_values.append(best_fitness)\n",
    "    return best_individual"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "As point of reference, we will compare against our initial default values again:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "fitness_values = []\n",
    "mutation_values = []\n",
    "max_evaluations = 10000\n",
    "\n",
    "configuration['P_xover'] = 0.7\n",
    "configuration['P_mutation'] = 0.07\n",
    "configuration['population_size'] = 100\n",
    "configuration['tournament_size'] = 2\n",
    "\n",
    "controlled_ga()\n",
    "controlled_fitness = fitness_values\n",
    "\n",
    "fitness_values = []\n",
    "ga()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "outputs": [],
   "source": [
    "plt.plot(mutation_values, label=\"Mutation Probability\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "plt.plot(controlled_fitness, label=\"Controlled\")\n",
    "plt.plot(fitness_values, label=\"Regular\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Adaptive Evolutionary Algorithms"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "The function that we used to control the mutation rate depends only on the phase of the search, but does not consider how the search is performing at any particular point in time. The idea of _adaptive_ parameter control is to adjust parameters based on how well the search is going. As a basic example, if mutation leads to more fitness improvement than crossover, we may want to increase the probability of doing mutation, and reduce the probability of doing crossover, and vice versa. In order to do this, we need to keep track of the fitness values before and after applying these operators:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def adaptive_crossover(parent1, parent2):\n",
    "    fitness_pre = get_fitness(parent1) + get_fitness(parent2)\n",
    "    offspring1, offspring2 = crossover(parent1, parent2)\n",
    "    fitness_post = get_fitness(offspring1) + get_fitness(offspring2)\n",
    "    xover_performance.append(fitness_post - fitness_pre)\n",
    "    return offspring1, offspring2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def adaptive_mutate(individual, P_mutate):\n",
    "    fitness_pre = get_fitness(individual)\n",
    "    mutated = mutate(individual, P_mutate)\n",
    "    fitness_post = get_fitness(mutated)\n",
    "    mutation_performance.append(fitness_post - fitness_pre)\n",
    "    return mutated"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Note that, technically, we are introducing additional fitness evaluations. However, to some degree these could be avoided (though not all) using some proper caching of fitness values."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "mutation_performance = []\n",
    "xover_performance = []\n",
    "\n",
    "adaptive_values = []\n",
    "xover_values = []\n",
    "mutation_values = []"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "The individual evolution step is the same as usual, we just need to make sure it uses the correct parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def adaptive_evolution_step(population, mutation_probability, xover_probability):\n",
    "    new_population = []\n",
    "    \n",
    "    while len(new_population) < len(population):\n",
    "        parent1 = selection(population)\n",
    "        parent2 = selection(population)\n",
    "\n",
    "        if random.random() < xover_probability:\n",
    "            offspring1, offspring2 = adaptive_crossover(parent1, parent2)\n",
    "        else:\n",
    "            offspring1, offspring2 = parent1, parent2\n",
    "\n",
    "        offspring1 = adaptive_mutate(offspring1, mutation_probability)\n",
    "        offspring2 = adaptive_mutate(offspring2, mutation_probability)\n",
    "\n",
    "        new_population.append(offspring1)\n",
    "        new_population.append(offspring2)\n",
    "\n",
    "    population.clear()\n",
    "    population.extend(new_population)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "After each step we need to re-evaluate the performance of the operators, and decide whether and how to update the parameters."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def adaptive_ga():\n",
    "    population = [get_random_solution() for _ in range(configuration[\"population_size\"])]\n",
    "    best_individual = max(population, key=lambda k: get_fitness(k))\n",
    "    best_fitness = get_fitness(best_individual)\n",
    "    print(f\"Initial population, best fitness: {best_fitness}\")\n",
    "\n",
    "    iteration = 0\n",
    "    max_iterations = (max_evaluations / configuration[\"population_size\"])\n",
    "    xover_probability = 0.7\n",
    "    mutation_probability = 0.07\n",
    "    delta = 0.01\n",
    "    \n",
    "    while iteration < max_iterations:\n",
    "        mutation_values.append(mutation_probability)\n",
    "        xover_values.append(xover_probability)\n",
    "        fitness_values.append(best_fitness)\n",
    "        iteration += 1\n",
    "        adaptive_evolution_step(population, mutation_probability, xover_probability)\n",
    "        \n",
    "        current_best = max(population, key=lambda k: get_fitness(k))\n",
    "        current_fitness = get_fitness(current_best)\n",
    "        if current_fitness > best_fitness:\n",
    "            print(f\"Iteration {iteration}, best fitness: {best_fitness}\")\n",
    "            best_individual = copy.deepcopy(current_best)\n",
    "            best_fitness = current_fitness\n",
    "            \n",
    "        max_mutation = max(mutation_performance) if mutation_performance else 0\n",
    "        max_xover = max(xover_performance) if xover_performance else 0\n",
    "\n",
    "        if max_mutation < max_xover:\n",
    "            xover_probability = min(xover_probability + delta, 1.0)\n",
    "            mutation_probability = max(mutation_probability - delta, 0.01)\n",
    "        elif max_xover < max_mutation:\n",
    "            xover_probability = max(xover_probability - delta, 0.0)\n",
    "            mutation_probability = min(mutation_probability + delta, 1.0)\n",
    "\n",
    "        mutation_performance.clear()\n",
    "        xover_performance.clear()\n",
    "\n",
    "    fitness_values.append(best_fitness)\n",
    "    return best_individual"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "We'll compare how this performs compared to our initial default configuration again."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "fitness_values = []\n",
    "num_mutations = []\n",
    "\n",
    "configuration['P_xover'] = 0.7\n",
    "configuration['P_mutation'] = 0.07\n",
    "configuration['population_size'] = 100\n",
    "configuration['tournament_size'] = 2\n",
    "\n",
    "mutation_values = []\n",
    "xover_values = []\n",
    "adaptive_ga()\n",
    "adaptive_fitness = fitness_values\n",
    "\n",
    "fitness_values = []\n",
    "ga()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "plt.plot(mutation_values, label=\"Mutation\")\n",
    "plt.plot(xover_values, label=\"Xover\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "plt.plot(adaptive_fitness, label=\"Adaptive\")\n",
    "plt.plot(fitness_values, label=\"Regular\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Self-Adaptive Evolutionary Algorithms"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Adaptive evolutionary algorithms control some parameters, but thereby introduce new ones (e.g., when and how and how much to update the controlled parameters). What should we set these new parameters to? And should these values be constant, or also adaptive? To avoid this problem altogether we can just let our algorithm decide on its own what values to use, and turn it into a _self-adaptive_ evolutionary algorithm."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "We thus add the parameters we want the algorithm to adapt to the representation:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def get_random_solution():\n",
    "    solution = L([])\n",
    "    for test in range(num_tests):\n",
    "        solution.append(get_random_test())\n",
    "\n",
    "    solution.probabilities = {\"P_xover\" : random.random(),\n",
    "                              \"P_mutation\" : random.random()}\n",
    "    return solution"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "When crossing over two individuals, for simplicity we will simply let them inherit their parents' parameter values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def self_crossover(parent1, parent2):\n",
    "    offspring1, offspring2 = crossover(parent1, parent2)\n",
    "\n",
    "    # Make sure copies have the probabilities\n",
    "    offspring1.probabilities = parent1.probabilities.copy()\n",
    "    offspring2.probabilities = parent2.probabilities.copy()\n",
    "\n",
    "    return offspring1, offspring2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Now comes the important bit: When mutating individuals, we not only alter the values (i.e., the test data in our search problem), but also the parameter values. We will add some Gaussian noise to mutate parameters by a little bit, so here's a helper function to keep values in bounds when applying the Gaussian noise:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def crop(value, boundary_min, boundary_max):\n",
    "    if value > boundary_max:\n",
    "        return boundary_max\n",
    "    elif value < boundary_min:\n",
    "        return boundary_min\n",
    "    else:\n",
    "        return value"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "The actual mutation just calls our usual mutation, and then applies the Gaussian noise -- using the current mutation probability to decide whether or not to mutate each parameter:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def self_mutate(individual, P_mutate):\n",
    "    \n",
    "    mutated = mutate(individual, P_mutate)\n",
    "    mutated.probabilities = individual.probabilities.copy()\n",
    "    \n",
    "    for p in mutated.probabilities.keys():\n",
    "        if random.random() < P_mutate:\n",
    "            mutated.probabilities[p] = crop(random.gauss(mutated.probabilities[p], 0.05), 0.001, 1)\n",
    "            \n",
    "    return mutated"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "When given two parents, we use the average of the crossover probabilities of the parents as our crossover probability; when mutating an individual, we use its own mutation probability:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def self_evolution_step(population):\n",
    "    new_population = []\n",
    "    while len(new_population) < len(population):\n",
    "        parent1 = selection(population)\n",
    "        parent2 = selection(population)\n",
    "\n",
    "        P_xover = (parent1.probabilities[\"P_xover\"] + parent2.probabilities[\"P_xover\"])/2\n",
    "        if random.random() < P_xover:\n",
    "            offspring1, offspring2 = self_crossover(parent1, parent2)\n",
    "        else:\n",
    "            offspring1, offspring2 = parent1, parent2\n",
    "\n",
    "        offspring1 = self_mutate(offspring1, offspring1.probabilities[\"P_mutation\"])\n",
    "        offspring2 = self_mutate(offspring2, offspring2.probabilities[\"P_mutation\"])\n",
    "\n",
    "        new_population.append(offspring1)\n",
    "        new_population.append(offspring2)\n",
    "\n",
    "    population.clear()\n",
    "    population.extend(new_population)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "The evolutionary algorithm itself is unchanged, but we add some code to keep track of the average parameters in the population so we can look at that after the evolution:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def self_adaptive_ga():\n",
    "    population = [get_random_solution() for _ in range(configuration[\"population_size\"])]\n",
    "    best_individual = max(population, key=lambda k: get_fitness(k))\n",
    "    best_fitness = get_fitness(best_individual)\n",
    "    print(f\"Initial population, best fitness: {best_fitness}\")\n",
    "\n",
    "    iteration = 0\n",
    "    max_iterations = (max_evaluations / configuration[\"population_size\"])\n",
    "    \n",
    "    while iteration < max_iterations:\n",
    "        fitness_values.append(best_fitness)\n",
    "        mutation_values.append(sum([p.probabilities[\"P_mutation\"] for p in population])/len(population))\n",
    "        xover_values.append(sum([p.probabilities[\"P_xover\"] for p in population])/len(population))\n",
    "        iteration += 1\n",
    "        self_evolution_step(population)\n",
    "        \n",
    "        current_best = max(population, key=lambda k: get_fitness(k))\n",
    "        current_fitness = get_fitness(current_best)\n",
    "        if current_fitness > best_fitness:\n",
    "            print(f\"Iteration {iteration}, best fitness: {best_fitness}\")\n",
    "            best_individual = copy.deepcopy(current_best)\n",
    "            best_fitness = current_fitness\n",
    "            \n",
    "    fitness_values.append(best_fitness)\n",
    "    return best_individual"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Once again we compare against our initial default values."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "max_evaluations = 10000\n",
    "fitness_values = []\n",
    "configuration['P_xover'] = 0.7\n",
    "configuration['P_mutation'] = 0.07\n",
    "configuration['population_size'] = 100\n",
    "configuration['tournament_size'] = 2\n",
    "\n",
    "self_values = []\n",
    "xover_values = []\n",
    "mutation_values = []\n",
    "\n",
    "self_adaptive_ga()\n",
    "self_values = fitness_values\n",
    "\n",
    "fitness_values = []\n",
    "ga()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "plt.plot(mutation_values, label=\"Mutation\")\n",
    "plt.plot(xover_values, label=\"Xover\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "plt.plot(self_values, label=\"Self-Adaptive\")\n",
    "plt.plot(fitness_values, label=\"Regular\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Generalisation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "Throughout this notebook, when we compared algorithms we only looked at individual runs; it's difficult to generalise from an individual run. Furthermore, we always applied the search to the exact same problem. In order to better compare the different configurations, we need to consider _different_ problems."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "from IPython.utils import io\n",
    "\n",
    "def run_times(algorithm, configuration, repetitions):\n",
    "    global fitness_values\n",
    "    result = []\n",
    "    for i in range(repetitions):\n",
    "        fitness_values = []\n",
    "        with io.capture_output() as captured: \n",
    "            algorithm()\n",
    "        result.append(fitness_values[-1])\n",
    "    return result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "parameters = { 0 : [\"a\", \"b\", \"c\", \"d\"],\n",
    "               1 : [\"d\", \"e\", \"f\", \"g\", \"x\"],\n",
    "               2 : [\"g\", \"h\", \"i\"]}\n",
    "num_tests = 20\n",
    "configuration['P_xover'] = 0.7\n",
    "configuration['P_mutation'] = 0.07\n",
    "configuration['population_size'] = 100\n",
    "configuration['tournament_size'] = 2\n",
    "results = {\n",
    "    \"Default\"       : run_times(ga, configuration, 10),\n",
    "    \"RSM\"           : run_times(ga, rsm_config, 10),\n",
    "    \"Controlled\"    : run_times(controlled_ga, configuration, 10),\n",
    "    \"Adaptive\"      : run_times(adaptive_ga, configuration, 10),\n",
    "    \"Self-Adaptive\" : run_times(self_adaptive_ga, configuration, 10)\n",
    "}\n",
    "fig, ax = plt.subplots()\n",
    "ax.boxplot(results.values())\n",
    "ax.set_xticklabels(results.keys())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "parameters = { 0 : [\"a\", \"b\", \"c\", \"d\", \"e\", \"f\", \"g\"],\n",
    "               1 : [\"d\", \"e\", \"f\", \"g\", \"x\"],\n",
    "               2 : [\"d\", \"e\", \"f\", \"g\"],\n",
    "               3 : [\"d\", \"e\", \"f\"],\n",
    "               4 : [\"g\", \"h\"]}\n",
    "num_tests = 30\n",
    "configuration['P_xover'] = 0.7\n",
    "configuration['P_mutation'] = 0.07\n",
    "configuration['population_size'] = 100\n",
    "configuration['tournament_size'] = 2\n",
    "results = {\n",
    "    \"Default\"       : run_times(ga, configuration, 10),\n",
    "    \"RSM\"           : run_times(ga, rsm_config, 10),\n",
    "    \"Controlled\"    : run_times(controlled_ga, configuration, 10),\n",
    "    \"Adaptive\"      : run_times(adaptive_ga, configuration, 10),\n",
    "    \"Self-Adaptive\" : run_times(self_adaptive_ga, configuration, 10)\n",
    "}\n",
    "fig, ax = plt.subplots()\n",
    "ax.boxplot(results.values())\n",
    "ax.set_xticklabels(results.keys())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "parameters = { 0 : [\"0\", \"1\"],\n",
    "               1 : [\"0\", \"1\"],\n",
    "               2 : [\"0\", \"1\"],\n",
    "               3 : [\"0\", \"1\"],\n",
    "               4 : [\"0\", \"1\"],\n",
    "               5 : [\"0\", \"1\"],\n",
    "               6 : [\"0\", \"1\"],\n",
    "               7 : [\"0\", \"1\"],\n",
    "               8 : [\"0\", \"1\"],\n",
    "               9 : [\"0\", \"1\"]}\n",
    "num_tests = 6 # 10 is required\n",
    "configuration['P_xover'] = 0.7\n",
    "configuration['P_mutation'] = 0.07\n",
    "configuration['population_size'] = 100\n",
    "configuration['tournament_size'] = 2\n",
    "results = {\n",
    "    \"Default\"       : run_times(ga, configuration, 10),\n",
    "    \"RSM\"           : run_times(ga, rsm_config, 10),\n",
    "    \"Controlled\"    : run_times(controlled_ga, configuration, 10),\n",
    "    \"Adaptive\"      : run_times(adaptive_ga, configuration, 10),\n",
    "    \"Self-Adaptive\" : run_times(self_adaptive_ga, configuration, 10)\n",
    "}\n",
    "fig, ax = plt.subplots()\n",
    "ax.boxplot(results.values())\n",
    "ax.set_xticklabels(results.keys())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "## Outlook: Hyper-Heuristics"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "In this chapter we considered how to tune and adjust the quantitative parameter values for a given search algorithm. However, sometimes different qualitative parameters (e.g., different types of mutation operators) may be beneficial at different times throughout the search, and sometimes even entirely different algorithms may perform better at different points during the search. Indeed there might not be a single optimal heuristic for a given problem type. The idea of _hyper heuristics_ is to apply different heuristics to different parts or phases of the solution process.\n",
    "A hyper-heuristic is a search method for selecting or generating heuristics to solve computational search problems, and so it can, essentially, produce entirely new and adaptive meta-heuristic search algorithms."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "An example for a hyper-heuristic application in software engineering can be found in the following article:\n",
    "\n",
    "- Jia, Y., Cohen, M. B., Harman, M., & Petke, J. Learning combinatorial interaction test generation strategies using hyperheuristic search. In Proceedings of the 37th International Conference on Software Engineering-Volume 1 (pp. 540-550). IEEE Press, 2015\n",
    "\n",
    "In this article, a selective hyper-heuristic is applied to combinatorial interaction testing, i.e. the example problem we used in this notebook. Simulated annealing is the central meta-heuristic, but six different neighbourhood exploration operators are defined, and reinforcement learning is applied in order to select which of them should be applied at each iteration of the simulated annealing algorithm."
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Slideshow",
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.5"
  },
  "rise": {
   "enable_chalkboard": true,
   "scroll": true
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}