Chapter05_MultiObjectiveOpt.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "e9ee845e-f52e-4382-bf6b-aa655d9b2132",
   "metadata": {
    "tags": []
   },
   "source": [
    "<h1>Multi-Objective Optimization</h1>\n",
    "It happens to have multiple objectives for the optimization problem. <br>\n",
    "<img src=\"images/moo_form.png\" width=300>\n",
    "Consider buying a laptop, perhaps you want it to be powerful and also you want to pay less money. So, you have two objectives which might contradict each others. In multi-objective problems, we are looking for the pareto-frontier of the solution. Any point on the pareto is not dominated by other solutions. Check the following figure for clarification<br>\n",
    "<img src=\"images/pareto_frontier.png\">\n",
    "<br>\n",
    "There are a variety of packages for multi-objective optimization in Python. Here, we introduce pymoo.<br>\n",
    "<img src=\"images/pymoo_logo.png\" width=300>\n",
    "Pymoo is an open source package for multi-objective optimization in python. It supports a variety of gradient-free optimization techniques.<br>\n",
    "For installation use: <code>conda install -c conda-forge pymoo</code><br>\n",
    "Let's review some of the capabilities here.<br><br>\n",
    "Pymoo.factory has various prebuilt items. For example, Pymoo has all benchmark functions predefined.\n",
    "<h2>Visualizing Benchmark Objective Functions</h2>\n",
    "Rosenberk function: Single Objective Benchmark function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6db001af-9cfa-4784-8ea1-0a8fe8f328aa",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from pymoo.factory import get_problem, get_visualization\n",
    "\n",
    "rosen_problem = get_problem(\"rosenbrock\", n_var=2)\n",
    "_=get_visualization(\"fitness-landscape\", rosen_problem, angle=(45, 45), _type=\"surface\").show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4e710e33-389a-4d24-b695-77f6d8c8cb60",
   "metadata": {},
   "source": [
    "<a href=\"https://www.researchgate.net/profile/Thanh-Binh-To/publication/2446876_MOBES_A_Multiobjective_Evolution_Strategy_for_Constrained_Optimization_Problems/links/53eb422e0cf28f342f45251b/MOBES-A-Multiobjective-Evolution-Strategy-for-Constrained-Optimization-Problems.pdf\"> \n",
    "Binh and Korn</a> problem is a benchmark for multi-objective optimization problem. It is defined as:<br>\n",
    "<img src=\"images/bnh.png\" width=300>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "25052d79-7b07-4d02-bf3e-ae47a98d5905",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pymoo.util.plotting import plot\n",
    "bnh_problem = get_problem(\"bnh\")\n",
    "plot(bnh_problem.pareto_front(), no_fill=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7848c0bc-f9f2-462e-a512-6d584e9a9492",
   "metadata": {},
   "source": [
    "<h2>Solving an optimization Problem (Single Objective)</h2>\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8d377ce1-fa32-481a-9166-aa63292afdfe",
   "metadata": {},
   "outputs": [],
   "source": [
    "# solving a problem: Here we solve a benchmark problem\n",
    "# Step 1: Load the solver and needed operators\n",
    "from pymoo.algorithms.soo.nonconvex.ga import GA\n",
    "from pymoo.factory import get_termination\n",
    "# define the solver\n",
    "algorithm = GA(pop_size=30)\n",
    "# Defining a stop condition\n",
    "termination = get_termination(\"n_gen\", 40)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f80f4407-e2bb-44e8-aeb5-1b1b57c10bac",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Now minimizing the objective\n",
    "from pymoo.optimize import minimize\n",
    "\n",
    "res = minimize(rosen_problem,\n",
    "               algorithm,\n",
    "               termination,\n",
    "               seed=1,)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8888ea10-d51e-4089-ba00-988a21831140",
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"X_opt:\", res.X, \", Function Value at x_opt:\", res.F)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7991b5b3-8cbc-47e5-b84a-afe9829f6d43",
   "metadata": {},
   "source": [
    "<img src=\"https://img.icons8.com/color/344/light.png\" width=80 alt=\"Hint\"> __hint:__ \n",
    "The result object has various useful fields:\n",
    "<ul class=\"simple\">\n",
    "<li><p><code class=\"docutils literal notranslate\"><span class=\"pre\">res.X</span></code>: Design space values are</p></li>\n",
    "<li><p><code class=\"docutils literal notranslate\"><span class=\"pre\">res.F</span></code>: Objective spaces values</p></li>\n",
    "<li><p><code class=\"docutils literal notranslate\"><span class=\"pre\">res.G</span></code>: Constraint values</p></li>\n",
    "<li><p><code class=\"docutils literal notranslate\"><span class=\"pre\">res.CV</span></code>: Aggregated constraint violation</p></li>\n",
    "<li><p><code class=\"docutils literal notranslate\"><span class=\"pre\">res.algorithm</span></code>: Algorithm object which has been iterated over</p></li>\n",
    "<li><p><code class=\"docutils literal notranslate\"><span class=\"pre\">res.opt</span></code>: The solutions as a <code class=\"docutils literal notranslate\"><span class=\"pre\">Population</span></code> object.</p></li>\n",
    "<li><p><code class=\"docutils literal notranslate\"><span class=\"pre\">res.pop</span></code>: The final Population</p></li>\n",
    "<li><p><code class=\"docutils literal notranslate\"><span class=\"pre\">res.history</span></code>: The history of the algorithm. (only if <code class=\"docutils literal notranslate\"><span class=\"pre\">save_history</span></code> has been enabled during the algorithm initialization)</p></li>\n",
    "<li><p><code class=\"docutils literal notranslate\"><span class=\"pre\">res.time</span></code>: The time required to run the algorithm</p></li>\n",
    "</ul>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "97dedc87-2159-4627-8353-ead5b2224dd1",
   "metadata": {},
   "source": [
    "<img src=\"https://img.icons8.com/external-soft-fill-juicy-fish/344/external-maths-school-soft-fill-soft-fill-juicy-fish.png\" alt=\"Math Tip\" width=50>\n",
    "<h5>Genetic Algorithm Review</h5>\n",
    "The solver that we used in above example is an advanced version of Gnetic Algorithm (GA). Here we review the basic GA algorithm for integer optimization. \n",
    "<table><tr><td><img src=\"images/genetic_algo.jpeg\" width=400></td><td>\n",
    "<ul>\n",
    "    <li>Each parameter is called a gene.</li>\n",
    "    <li>The sequence of all parameters is called chromosome. </li>\n",
    "    <li> The objective function at each chromosome has a value</li>\n",
    "    <li> In each iteration is called a generation.</li>\n",
    "    <li>In each generation, there are a population of chromosome.<li>\n",
    "    <li>The population is initialized randomly.</li>\n",
    "    <li>We pass the best chromosomes to the next generation</li>\n",
    "    <li> In each generation there are new chromosomes obtained by cross-over and mutation operations</li>\n",
    "    <li>Cross-over has two parents and they are merged to make the offsprings</li>\n",
    "    <li>Mutation is selected randomly and one or more genes are updated</li>\n",
    "</ul>\n",
    "</td></tr></table>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c5ac7b84-96a0-4cd8-8075-11c1c6b01b74",
   "metadata": {},
   "source": [
    "<img src= \"https://img.icons8.com/external-flaticons-flat-flat-icons/344/external-question-100-most-used-icons-flaticons-flat-flat-icons.png\" width=70> __Question (Rastrigin Problem)__:<br>\n",
    "The following problem is called Rastrigin and is one optimization benchmarks.<br>\n",
    "<img src=\"images/rastrigin.png\" width=300><br>\n",
    "It has multiple local minimas and a global minima at (0,0). This problem is already implemented by factory. Please, write a program to do the followings:<br>\n",
    "<ul>\n",
    "    <li> Load the problem for <code>name=\"rastrigin\"</code> and <code>n_vars=2</code>.</li>\n",
    "    <li> Plot the contours of the objective function (use <code>_type=\"contour\"</code>).</li>\n",
    "    <li> The stop condition is iterating for 100 times. <code> n_gen=100</code></li>\n",
    "    <li> Solve the optimization problem using GA with <code>pop_size=300</code></li>\n",
    "    <li> Print the optimal point found and the value of the function at the optimal point. </li>\n",
    "</ul>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "33d1887f-962d-4a04-977b-4d8ec1afe058",
   "metadata": {},
   "outputs": [],
   "source": [
    "# let's run the last optimization with more verobosity\n",
    "res = minimize(rosen_problem,\n",
    "               algorithm,\n",
    "               termination,\n",
    "               seed=1,\n",
    "              verbose=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1ac22245-b106-4eb0-8438-bc5bcb61ca73",
   "metadata": {},
   "source": [
    "<img src=\"https://img.icons8.com/color/344/light.png\" width=80 alt=\"Hint\"> __hint:__ \n",
    "The display information contains useful diagnosis data and can provide intuition about what is going on. For checking this output you should set the verbose=True option.\n",
    "<table class=\"colwidths-given table table-responsive w-100 d-block d-md-table pure-table pure-table-bordered my-table\" id=\"id2\">\n",
    "<caption><span class=\"caption-text\">Types of Output</span><a class=\"headerlink\" href=\"#id2\" title=\"Permalink to this table\">¶</a></caption>\n",
    "<colgroup>\n",
    "<col style=\"width: 14%\">\n",
    "<col style=\"width: 86%\">\n",
    "</colgroup>\n",
    "<thead>\n",
    "<tr class=\"row-odd\"><th class=\"head\"><p>Name</p></th>\n",
    "<th class=\"head\"><p>Description</p></th>\n",
    "</tr>\n",
    "</thead>\n",
    "<tbody>\n",
    "<tr class=\"row-even\"><td><p><strong>n_gen</strong></p></td>\n",
    "<td><p>The current number of generations or iterations until this point.</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><p><strong>n_eval</strong></p></td>\n",
    "<td><p>The number of function evaluations so far.</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><p><strong>n_nds</strong></p></td>\n",
    "<td><p>For multi-objective problems, the number of <em>non-dominated</em> solutions of the optima found.</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><p><strong>cv (min)</strong></p></td>\n",
    "<td><p>The minimum constraint violation (CV) in the current population</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><p><strong>cv (avg)</strong></p></td>\n",
    "<td><p>The average constraint violation (CV) of the current population</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><p><strong>f_opt</strong></p></td>\n",
    "<td><p>For single-objective problems, the best function value found so far.</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><p><strong>f_gap</strong></p></td>\n",
    "<td><p>For single-objective problems, the best gap to the optimum (only printed if the optimum is <em>known</em>).</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><p><strong>eps/indicator</strong></p></td>\n",
    "<td><p>For multi-objective problems, the change of the indicator (ideal, nadir, f) over the last few generations (only printed if the Pareto-front is <em>unknown</em>). For more information we encourage you to have a look at the corresponding publication (<a class=\"bibtex reference internal\" href=\"../references.html#running\" id=\"id1\">[22]</a>, <a class=\"reference external\" href=\"https://www.egr.msu.edu/~kdeb/papers/c2020003.pdf\">pdf</a>).</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><p><strong>igd/gd/hv</strong></p></td>\n",
    "<td><p>For multi-objective problems, the performance indicator (only printed if the Pareto-front is <em>known</em>).</p></td>\n",
    "</tr>\n",
    "</tbody>\n",
    "</table>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c165d573-b0df-41cc-95f4-8d1bd70dd002",
   "metadata": {},
   "source": [
    "<h2> Solving Multi-objective Optimization Problem</h2>\n",
    "Let's solve BNH Problem"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "64c04d03-cc8c-4f66-a3df-032c683fb530",
   "metadata": {},
   "outputs": [],
   "source": [
    "# We need a multi-objective solver here\n",
    "from pymoo.algorithms.moo.nsga2 import NSGA2\n",
    "from pymoo.factory import get_sampling, get_crossover, get_mutation\n",
    "moo_algorithm = NSGA2(\n",
    "    pop_size=40,\n",
    "    n_offsprings=10,\n",
    "    sampling=get_sampling(\"real_random\"),\n",
    "    crossover=get_crossover(\"real_sbx\", prob=0.9, eta=15),\n",
    "    mutation=get_mutation(\"real_pm\", eta=20),\n",
    "    eliminate_duplicates=True\n",
    ")\n",
    "\n",
    "res = minimize(bnh_problem,\n",
    "               moo_algorithm,\n",
    "               termination,\n",
    "               seed=1,\n",
    "               save_history=True,\n",
    "               verbose=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "38be6b8f-4057-4b08-b107-5d40ff4e1785",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Let's plot pareto of the solution\n",
    "%matplotlib inline\n",
    "import matplotlib as mpl\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(10, 8))\n",
    "ax.plot(res.F[:,0], res.F[:,1], linestyle=\"\", marker=\".\", label='Solutions')  \n",
    "ax.set_xlabel('F_1 Objective') \n",
    "ax.set_ylabel('F_2 Objective')\n",
    "ax.set_title(\"Pareto Frontier\") \n",
    "ax.legend();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "be512b3c-e5e2-4509-93ee-935feca30bfe",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Let's say we check solutions and we are interested in the one which has the minimum |F1|+ |F2|\n",
    "import numpy as np\n",
    "selection = np.argmin(np.abs(res.F[:,0]) + np.abs(res.F[:,1]))\n",
    "print(\"Final Selected Solution: \\n\\tSelction Index: \",selection, \"X_opt: \",\n",
    "      res.X[selection,:], \", Objectives: F1->\", res.F[selection, 0], \", F2->\", res.F[selection, 1])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e9abff12-26bc-4d85-957c-7a4e26d5f122",
   "metadata": {},
   "source": [
    "<img src=\"https://img.icons8.com/external-soft-fill-juicy-fish/344/external-maths-school-soft-fill-soft-fill-juicy-fish.png\" alt=\"Math Tip\" width=50><font size=5>Ideal and Nadir Normalization</font><br>What we have done above is called __Multi-Criteria Decision Making (MCDM)__ in the optimization literature. It means, once you have the parto space, you can have another criteria to narrow down your choice of optimal point. <br>\n",
    "Beware of the fact that objectives could have different scale and you might want to rescale them before applying any logic. One common way of normalization is using ideal and nadir points. Ideal is a vector with minimum of each objective among all pareto points and nadir is the maximum. Once ideal and nadir points are calculated (they are estimation for the real points of ideal and nadir) one can normalize using ideal and nadir as below:<br>\n",
    "<img src=\"images/ideal_nadir.png\" width=400><br>\n",
    "Then a use might have some wishes or weights for each objective function. Then, one can use __Augmented Scalarization Function (ASF)__ metric to select a solution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2f82069b-a0ec-487d-9035-e218f978b097",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Let's do ideal nadir normalization this time\n",
    "approx_ideal = res.F.min(axis=0)\n",
    "approx_nadir = res.F.max(axis=0)\n",
    "# Normalized Pareto\n",
    "nF = (res.F - approx_ideal) / (approx_nadir - approx_ideal)\n",
    "# plotting pareto again and show ideal and nadir\n",
    "fig, ax = plt.subplots(figsize=(10, 8))\n",
    "ax.plot(nF[:,0], nF[:,1], linestyle=\"\", marker=\".\", label='Normalized Solutions')  \n",
    "ax.set_xlabel('Normalized F_1 Objective') \n",
    "ax.set_ylabel('Normalized F_2 Objective')\n",
    "ax.set_title(\"Normalized Pareto Frontier\") \n",
    "ax.legend();\n",
    "print(\"Ideal: \", approx_ideal)\n",
    "print(\"Nadir: \", approx_nadir)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fb36008b-d3a1-4ff9-ad6b-4473e58f1349",
   "metadata": {},
   "outputs": [],
   "source": [
    "weights = np.array([0.2, 0.8])\n",
    "from pymoo.decomposition.asf import ASF\n",
    "\n",
    "decomp = ASF()\n",
    "selection = decomp.do(nF, 1/weights).argmin()\n",
    "print(\"Final Selected Solution: \\n\\tSelction Index: \",selection, \"X_opt: \",\n",
    "      res.X[selection,:], \", Objectives: F1->\", res.F[selection, 0], \", F2->\", res.F[selection, 1])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2ca5ef53-8554-4328-a463-1aef76623175",
   "metadata": {},
   "source": [
    "<img src= \"https://img.icons8.com/external-flaticons-flat-flat-icons/344/external-question-100-most-used-icons-flaticons-flat-flat-icons.png\" width=70> __Question (Tanaka Problem)__:<br>\n",
    "The following multi-objective problem is called Tanaka problem. It is pre-built by the factory module <code> name =\"tnk\"</code>. <br>\n",
    "<img src=\"images/tanaka.png\" width=400><br>\n",
    "<ul>\n",
    "    <li>Load the program from the factory.</li>\n",
    "    <li>Display the answer (optimal pareto) without solving it.</li>\n",
    "    <li>Solve the problem using NSGAII with the following settings:\n",
    "        <code>pop_size=60, n_offsprings=10, crossover=get_crossover(\"real_sbx\", prob=0.7, eta=15), mutation=get_mutation(\"real_pm\", eta=20)</code></li>\n",
    "    <li>Plot the pareto of your solution</li>\n",
    "    <li>Apply ASF with weights of <code>[0.7, 0.3]</code> for MCDM</li>\n",
    "    <li>Print the outcome</li>\n",
    "    </ul>\n",
    "        "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "add12f79-2d30-4c52-a2e1-a4605c56927a",
   "metadata": {
    "tags": []
   },
   "source": [
    "<h2>Particle Swarm Optimization</h2>\n",
    "So far we used GA for both SOO and MOO. But pymoo covers a variety of algorithms. <code>pymoo.algorithms.soo</code> contains a variety of single objective algorithms and <code>pymoo.algorithms.moo</code> contais various multi objective algorithms. PSO is inspired by the way birds fly in groups. Each particle is associated with a velocity vector and position which are updated at each iteration. For the next positon, the particle move along its velocity vector. The velcity itself has three main components. <br>\n",
    "<table><tr><td><img src=\"images/pso.png\" width=350></td><td><img src=\"images/pso_update.png\" width=450></td><td><img src=\"https://upload.wikimedia.org/wikipedia/commons/e/ec/ParticleSwarmArrowsAnimation.gif\" width=400></td></tr></table>\n",
    "Check this <a href=\"http://pymoo.org/animations/pso.mp4\"> video</a> for an animation of the method.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9c29b478-cdf8-4e45-a273-aa5015d11323",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pymoo.algorithms.soo.nonconvex.pso import PSO\n",
    "\n",
    "problem = get_problem(\"ackley\")\n",
    "\n",
    "algorithm = PSO(max_velocity_rate=0.025, pop_size=100)\n",
    "\n",
    "res = minimize(problem,\n",
    "               algorithm,\n",
    "               seed=1,\n",
    "               save_history=True,\n",
    "               verbose=False)\n",
    "\n",
    "print(\"Best solution found: \\nX = %s\\nF = %s\" % (res.X, res.F))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c46ff6dd-1327-4186-b254-ac6024b75354",
   "metadata": {
    "tags": []
   },
   "source": [
    "<h3> Defining Custom Optimization Problem (SOO)</h3>\n",
    "So far we used benchmark problems. If you want to implement your own problem, you have to define it as class. Below, we implement Spehere problem (a benchmark problem) ourself."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4bfa8b0c-98a5-4f45-914e-49b79f608e8c",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pymoo.core.problem import Problem # The base class for problems\n",
    "\n",
    "class SphereWithConstraint(Problem):\n",
    "\n",
    "    def __init__(self):\n",
    "        super().__init__(n_var=10, n_obj=1, n_constr=1, xl=0.0, xu=1.0)\n",
    "\n",
    "    def _evaluate(self, x, out, *args, **kwargs):\n",
    "        out[\"F\"] = np.sum((x - 0.5) ** 2, axis=1)\n",
    "        out[\"G\"] = 0.1 - out[\"F\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "355b60c4-0804-4a5f-92f6-dc6420ac78f5",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Now let's solve it using PSO\n",
    "problem = SphereWithConstraint()\n",
    "res = minimize(problem,\n",
    "               algorithm,\n",
    "               termination,\n",
    "               seed=1,\n",
    "               save_history=True,\n",
    "               verbose=False)\n",
    "print(\"Best solution found: \\nX = %s\\nF = %s\" % (res.X, res.F))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e31d004a-0b40-428a-afce-7c265bb1580b",
   "metadata": {},
   "source": [
    "<table class=\"colwidths-given table table-responsive w-100 d-block d-md-table pure-table pure-table-bordered my-table\" id=\"id1\">\n",
    "<caption><span class=\"caption-text\">Problem Definition</span><a class=\"headerlink\" href=\"#id1\" title=\"Permalink to this table\">¶</a></caption>\n",
    "<colgroup>\n",
    "<col style=\"width: 14%\">\n",
    "<col style=\"width: 86%\">\n",
    "</colgroup>\n",
    "<thead>\n",
    "<tr class=\"row-odd\"><th class=\"head\"><p>Argument</p></th>\n",
    "<th class=\"head\"><p>Description</p></th>\n",
    "</tr>\n",
    "</thead>\n",
    "<tbody>\n",
    "<tr class=\"row-even\"><td><p><code class=\"code docutils literal notranslate\"><span class=\"pre\">n_var</span></code></p></td>\n",
    "<td><p>Integer value representing the number of design variables.</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><p><code class=\"code docutils literal notranslate\"><span class=\"pre\">n_obj</span></code></p></td>\n",
    "<td><p>Integer value representing the number of objectives.</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><p><code class=\"code docutils literal notranslate\"><span class=\"pre\">n_constr</span></code></p></td>\n",
    "<td><p>Integer value representing the number of constraints.</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><p><code class=\"code docutils literal notranslate\"><span class=\"pre\">xl</span></code></p></td>\n",
    "<td><p>Float or <code class=\"code docutils literal notranslate\"><span class=\"pre\">np.ndarray</span></code> of length <code class=\"code docutils literal notranslate\"><span class=\"pre\">n_var</span></code> representing the lower bounds of the design variables.</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-even\"><td><p><code class=\"code docutils literal notranslate\"><span class=\"pre\">xu</span></code></p></td>\n",
    "<td><p>Float or <code class=\"code docutils literal notranslate\"><span class=\"pre\">np.ndarray</span></code> of length <code class=\"code docutils literal notranslate\"><span class=\"pre\">n_var</span></code> representing the upper bounds of the design variables.</p></td>\n",
    "</tr>\n",
    "<tr class=\"row-odd\"><td><p><code class=\"code docutils literal notranslate\"><span class=\"pre\">type_var</span></code></p></td>\n",
    "<td><p>(optional) A type hint for the user what variable should be optimized.</p></td>\n",
    "</tr>\n",
    "</tbody>\n",
    "</table>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c28575ff-e3be-4415-be76-0b7be1aa96ad",
   "metadata": {},
   "source": [
    "<h3> Solving Custom Optimization Problem (MOO)</h3>\n",
    "For multi-objective problem you need to instantiate your class from another class called <code>ElementwiseProblem</code>. Let's implement the following example.\n",
    "<img src=\"images/custom_problem.png\" width=400>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0cc6828d-b5be-41a4-903f-c766b0c6d1dc",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Step 1: Defining the problem\n",
    "import numpy as np\n",
    "from pymoo.core.problem import ElementwiseProblem # We can define elementwise objectives.\n",
    "\n",
    "# We should inehrit from ElementwiseProblem\n",
    "class MyProblem(ElementwiseProblem):\n",
    "    def __init__(self):\n",
    "        super().__init__(n_var=2, # Number of optimization variables\n",
    "                         n_obj=2, # Number of objectives\n",
    "                         n_constr=2, # Number of constraints\n",
    "                         xl=np.array([-2,-2]), # Lower bounds for variables\n",
    "                         xu=np.array([2,2])) # Upper bounds for variables\n",
    "\n",
    "    # In addition to above properties, it has to have the following function.\n",
    "    def _evaluate(self, x, out, *args, **kwargs):\n",
    "        # x is the varaible for calculation of objective\n",
    "        # out keeps the output value, it is a dict object.\n",
    "        f1 = 100 * (x[0]**2 + x[1]**2)\n",
    "        f2 = (x[0]-1)**2 + x[1]**2\n",
    "\n",
    "        g1 = 2*(x[0]-0.1) * (x[0]-0.9) / 0.18\n",
    "        g2 = - 20*(x[0]-0.4) * (x[0]-0.6) / 4.8\n",
    "\n",
    "        out[\"F\"] = [f1, f2] # out[\"F\"] contains n_obj elements corresponding to objective values\n",
    "        out[\"G\"] = [g1, g2] # out[\"G\"] contains n_constr elements corresponding to contraint values of (G)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e002dfe1-e9ca-4d98-937e-1d9541a4ac26",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Step 2: Init an object and solve it\n",
    "problem = MyProblem()\n",
    "# Let's use the same solver to solve the problem\n",
    "res = minimize(problem,\n",
    "               moo_algorithm,\n",
    "               termination,\n",
    "               seed=1,\n",
    "               save_history=True,\n",
    "               verbose=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a41a4b9f-f34d-4bb6-b185-531fc98f69f7",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(10, 8))\n",
    "ax.plot(res.F[:,0], res.F[:,1], linestyle=\"\", marker=\".\", label='Solutions')  \n",
    "ax.set_xlabel('F_1 Objective') \n",
    "ax.set_ylabel('F_2 Objective')\n",
    "ax.set_title(\"Pareto Frontier\") \n",
    "ax.legend();"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a0d4bdef-5e3a-415a-938d-5569f757e66e",
   "metadata": {},
   "source": [
    "<img src= \"https://img.icons8.com/external-flaticons-flat-flat-icons/344/external-question-100-most-used-icons-flaticons-flat-flat-icons.png\" width=70> __Question (Customized Optimization)__:<br>\n",
    "This <a href=\"https://en.wikipedia.org/wiki/Test_functions_for_optimization\">wikipedia page</a> contains a list of benchmark functions for optimization. Look at this page and find <em>Goldstein-price</em> problem. <br>\n",
    "<ul>\n",
    "    <li> Implement a class for this problem. make sure all bounds are right</li>\n",
    "    <li> Solve the problem using PSO</li>\n",
    "    <li> Solve the problem using GA</li>\n",
    "    <li> Print results.</li>\n",
    " </ul>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "88ec6515-e611-4207-99cf-a2fced644d6b",
   "metadata": {},
   "source": [
    "<h2>Many Objective Optimization</h2>\n",
    "A <em>Many Objective Optimization</em> problem is a multi-objective problem with many (more than 4) objectives. The complexity of problem becomes intractable in many objective problems and many algorithms lose their efficiency. Pymoo offers some of the state of art algorithms for addressing such scenarios.<br>\n",
    "The following optimization problem proposed by Deb, Thiele, Laumanns and Zitzler usually referred to as DTLZ1 problem:<br>\n",
    "<img src=\"images/DTLZ1.png\" width=350><br>\n",
    "The last 𝑘=(𝑛−𝑀+1) variables are represented as 𝐱𝑀. Also g(.) is defined as:\n",
    "<img src=\"images/DTLZ1_cond.png\" width=350>\n",
    "This problem is considered as a benchmark for many objective optimization problem. Pymoo has this benchmark problem. Let's have a look at it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fd1e7b6e-cbf8-485c-9b88-ac6c6007a53e",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pymoo.factory import  get_reference_directions\n",
    "\n",
    "ref_dirs = get_reference_directions(\"das-dennis\", 3, n_partitions=12)\n",
    "\n",
    "problem = get_problem(\"dtlz1\")\n",
    "pf= problem.pareto_front(ref_dirs)\n",
    "_ = get_visualization(\"scatter\", angle=(45,45)).add(pf).show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cdd0b6aa-2ca0-42dd-bc42-e025c6d7ef12",
   "metadata": {},
   "source": [
    "<img src=\"https://img.icons8.com/external-soft-fill-juicy-fish/344/external-maths-school-soft-fill-soft-fill-juicy-fish.png\" alt=\"Math Tip\" width=50> <font size=4>Reference Directions:</font> \n",
    "Most evolutionary many-objective optimization (EMaO) algorithms, for instance NSGA3 or MOEAD, start with a description of a number of predefined set of reference directions on a unit simplex. Reference directions in an 𝑀-dimensional space are subject to<br>\n",
    "<center><img src=\"images/refdir_cond.png\" width=200></center>\n",
    "Das and Dennis’s structured approach for generating well-spaced reference points is a popular approach. The number of points on the unit simplex is determined by a parameter P (n_partitions), which indicates the number of gaps between two consecutive points along an objective axis. It turns out that the total number of points (n) on the unit simplex is: C(n+p-1, p). For better understanding the math of reference points, the reader is encouraged to have a look at the following paper. Here, we skip the details for sake of simplicity.\n",
    "\n",
    "\n",
    "Deb, Kalyanmoy, and J. Sundar. \"<a href=\"https://www.egr.msu.edu/~kdeb/papers/k2005012.pdf\">Reference point based multi-objective optimization using evolutionary algorithms.</a>\" Proceedings of the 8th annual conference on Genetic and evolutionary computation. 2006.\n",
    "\n",
    "<br>Now let's optimize the many objective problem. We use MOEAD (multi-objective evolutionary algorithm based on decomposition) Algorithm\n",
    "This algorithm is introduced in the following paper:\n",
    "Qingfu Zhang and Hui Li. <a href=\"https://link.springer.com/article/10.1007/s00500-014-1239-3#:~:text=normal%20boundary%20intersection-,The%20Multiobjective%20Evolutionary%20Algorithm%20based%20on%20Decomposition%20(MOEA%2FD),information%20from%20the%20optimization%20of\">A multi-objective evolutionary algorithm based on decomposition</a>. IEEE Transactions on Evolutionary Computation, Accepted, 2007.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c7d875ce-8426-4448-8e8e-8d2e8d67e625",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pymoo.algorithms.moo.moead import MOEAD\n",
    "algorithm = MOEAD(\n",
    "    ref_dirs,\n",
    "    n_neighbors=15,\n",
    "    prob_neighbor_mating=0.7,\n",
    ")\n",
    "\n",
    "res = minimize(problem,\n",
    "               algorithm,\n",
    "               ('n_gen', 200),\n",
    "               seed=1,\n",
    "               verbose=False)\n",
    "\n",
    "_ = get_visualization(\"scatter\").add(res.F).show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1e509b68-097c-44ac-9a6a-f11e2a550833",
   "metadata": {},
   "source": [
    "<img src=\"https://img.icons8.com/external-soft-fill-juicy-fish/344/external-maths-school-soft-fill-soft-fill-juicy-fish.png\" alt=\"Math Tip\" width=50><font size=4>MOEAD Algorithm</font><br> MOEA/D decompose the multi-objective optimization to many single objective optimization. It combines all objective functions using a linear combination to form a new single objective function. The minimum of this objective will be a point on the pareto front of the multi-objective problem. It repeats it for various combinations of weights to find the parte front.\n",
    "<img src=\"images/moead.png\" widht =400>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c0404f20-9bbc-45f1-a47a-d6278df20b7a",
   "metadata": {},
   "source": [
    "<img src=\"https://img.icons8.com/external-itim2101-lineal-color-itim2101/344/external-professor-life-style-avatar-itim2101-lineal-color-itim2101.png\" alt=\"Instructor\" width=60><br> Pymoo supports many single, multi and many objective algorithms. If you are interested check the following modules:<br>\n",
    "Single objective algorithms: <code>pymoo.algorithms.soo</code><br>\n",
    "For multi/many objective algorithms\" <code> pymoo.algorithms.moo</code>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "82909181-b16d-4b3c-869a-ae3968594fe9",
   "metadata": {},
   "source": [
    "<img src= \"https://img.icons8.com/external-flaticons-flat-flat-icons/344/external-question-100-most-used-icons-flaticons-flat-flat-icons.png\" width=70> __Question (Many Objective Optimization)__:<br>\n",
    "Write a program to do the followings:<br>\n",
    "<ul>\n",
    "    <li> Load the optimization benchmark problem called C1DTLZ1 (<code>name=\"c1dtlz1\"</code>)</li>\n",
    "    <li> Build the ref_dir as follow\n",
    "        <code>ref_dirs = get_reference_directions(\"das-dennis\", 3, n_partitions=12)</code></li>\n",
    "    <li> Before solving, plot the known pareto of the benchmark function</li>\n",
    "    <li> We are interested in solving the problem using an algorithm called C-TAEA (Two-Archive Evolutionary Algorithm for Constrained Multiobjective Optimization).For building an instance of the solver algorithm use the following snippet:<br>\n",
    "<code>from pymoo.algorithms.moo.ctaea import CTAEA\n",
    " algorithm = CTAEA(ref_dirs=ref_dirs)</code></li>\n",
    "    <li>Run the aglorithm for 600 generations pass <code>('n_gen', 600)</code> to minimize.</li>\n",
    "    <li>Time profile the algoritm using %timeit.</li>\n",
    "    </ul>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "25473755-4537-4b00-ab4e-12ca1e4a9583",
   "metadata": {},
   "source": [
    "<h2>Parallelization and Computation Performance</h2>\n",
    "Many of algorithms that we discussed are good for parallel computing. Each chromose in GA or each particle in PSO can be evaluated separately. It is important to implement the problem in the most efficient way. It is recommended to use numpy for writing the program and use <code>ElementwiseProblem</code> to write your problem. This improve the performance up to some degree. For more improvment you can use multi-threading or multiprocessing to run evaluations in parallel. \n",
    "Let's check an example"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5c8271b4-657a-4816-be15-348d429eb98a",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pymoo.core.problem import starmap_parallelized_eval\n",
    "from multiprocessing.pool import ThreadPool\n",
    "\n",
    "\n",
    "\n",
    "class MyProblem(ElementwiseProblem):\n",
    "\n",
    "    def __init__(self, **kwargs):\n",
    "        super().__init__(n_var=10, n_obj=1, n_constr=0, xl=-5, xu=5, **kwargs)\n",
    "\n",
    "    def _evaluate(self, x, out, *args, **kwargs):\n",
    "         out[\"F\"] = (x ** 8).sum() + (x ** 6).sum() + (x ** 4).sum() + (x ** 2).sum()\n",
    "\n",
    "# the number of threads to be used\n",
    "n_threads = 8\n",
    "# initialize the pool\n",
    "pool = ThreadPool(n_threads)\n",
    "\n",
    "# define two versions of the problem: paralellized and non-paralellized\n",
    "problem_non_par = MyProblem() \n",
    "problem_par = MyProblem(runner=pool.starmap, func_eval=starmap_parallelized_eval)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "310cedfe-7acb-4704-8570-f30b268aa070",
   "metadata": {},
   "outputs": [],
   "source": [
    "# solve it wiouth prallelization\n",
    "%timeit res = minimize(problem_non_par, PSO(pop_size=1000), seed=1, n_gen=100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5f1ec19a-63ef-479f-998e-02956c494ae6",
   "metadata": {},
   "outputs": [],
   "source": [
    "%timeit res = minimize(problem_par, PSO(pop_size=1000), seed=1, n_gen=100)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "870f9b8d-47e9-4a0b-9d79-f8919e72ded0",
   "metadata": {},
   "source": [
    "<img src=\"https://img.icons8.com/external-flaticons-lineal-color-flat-icons/344/external-coffee-cup-bakery-flaticons-lineal-color-flat-icons.png\" alt=\"icon\" width=80> <font size=5> Takehome: Callbacks</font>\n",
    "You can monitor the optimization problem by defining a callback function. For example it can be used to extract the learning curves."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9b8e9cd0-4a36-46e7-91b0-b630605d1116",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pymoo.core.callback import Callback\n",
    "\n",
    "class MyCallback(Callback):\n",
    "\n",
    "    def __init__(self) -> None:\n",
    "        super().__init__()\n",
    "        self.data[\"best\"] = []\n",
    "\n",
    "    def notify(self, algorithm):\n",
    "        self.data[\"best\"].append(algorithm.pop.get(\"F\").min())\n",
    "\n",
    "problem = get_problem(\"sphere\")\n",
    "algorithm = GA(pop_size=100)\n",
    "res = minimize(problem,\n",
    "               algorithm,\n",
    "               ('n_gen', 40),\n",
    "               seed=1,\n",
    "               callback=MyCallback(),\n",
    "               verbose=False)\n",
    "\n",
    "val = res.algorithm.callback.data[\"best\"]\n",
    "fig, ax = plt.subplots(figsize=(10, 8))\n",
    "ax.plot(np.arange(len(val)), val)\n",
    "ax.set_xlabel('Generation') \n",
    "ax.set_ylabel('Objective')\n",
    "_=ax.set_title(\"Learning Curve\") "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "76c8a2af-ffbd-43a7-a71e-15538ae7017b",
   "metadata": {},
   "source": [
    "<img src=\"https://img.icons8.com/color/344/light.png\" width=80 alt=\"Hint\"> __hint:__ We used the above example for demonstration. It can be done without callback using history:<br>\n",
    "<code>val = [e.opt.get(\"F\")[0] for e in res.history]</code>"
   ]
  }
 ],
 "metadata": {
  "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.9.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}