.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_numerical_methods/optimization/plot_ego.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note Click :ref:`here ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_numerical_methods_optimization_plot_ego.py: EfficientGlobalOptimization examples ==================================== .. GENERATED FROM PYTHON SOURCE LINES 6-18 The EGO algorithm (Jones, 1998) is an adaptative optimization method based on kriging. An initial design of experiment is used to build a first metamodel. At each iteration a new point that maximizes a criterion is chosen as optimizer candidate. The criterion uses a tradeoff between the metamodel value and the conditional variance. Then the new point is evaluated using the original model and the metamodel is relearnt on the extended design of experiment. .. GENERATED FROM PYTHON SOURCE LINES 21-30 .. code-block:: default import openturns as ot import openturns.viewer as viewer from matplotlib import pylab as plt import math as m ot.RandomGenerator.SetSeed(0) ot.ResourceMap.SetAsString("KrigingAlgorithm-LinearAlgebra", "LAPACK") ot.Log.Show(ot.Log.NONE) .. GENERATED FROM PYTHON SOURCE LINES 31-35 Ackley test-case ---------------- We first apply the EGO algorithm on the :ref:`Ackley function`. .. GENERATED FROM PYTHON SOURCE LINES 37-39 Define the problem ^^^^^^^^^^^^^^^^^^ .. GENERATED FROM PYTHON SOURCE LINES 41-42 The Ackley model is defined in the usecases module in a data class `AckleyModel` : .. GENERATED FROM PYTHON SOURCE LINES 42-45 .. code-block:: default from openturns.usecases import ackley_function as ackley_function am = ackley_function.AckleyModel() .. GENERATED FROM PYTHON SOURCE LINES 46-47 We get the Ackley function : .. GENERATED FROM PYTHON SOURCE LINES 47-49 .. code-block:: default model = am.model .. GENERATED FROM PYTHON SOURCE LINES 50-51 We specify the domain of the model : .. GENERATED FROM PYTHON SOURCE LINES 51-55 .. code-block:: default dim = am.dim lowerbound = [-4.0] * dim upperbound = [4.0] * dim .. GENERATED FROM PYTHON SOURCE LINES 56-57 We know that the global minimum is at the center of the domain. It is stored in the `AckleyModel` data class. .. GENERATED FROM PYTHON SOURCE LINES 57-59 .. code-block:: default xexact = am.x0 .. GENERATED FROM PYTHON SOURCE LINES 60-61 The minimum value attained `fexact` is : .. GENERATED FROM PYTHON SOURCE LINES 61-64 .. code-block:: default fexact = model(xexact) fexact .. raw:: html

[4.44089e-16]



.. GENERATED FROM PYTHON SOURCE LINES 65-69 .. code-block:: default graph = model.draw(lowerbound, upperbound, [100]*dim) graph.setTitle("Ackley function") view = viewer.View(graph) .. image:: /auto_numerical_methods/optimization/images/sphx_glr_plot_ego_001.png :alt: Ackley function :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 70-71 We see that the Ackley function has many local minimas. The global minimum, however, is unique and located at the center of the domain. .. GENERATED FROM PYTHON SOURCE LINES 73-77 Create the initial kriging ^^^^^^^^^^^^^^^^^^^^^^^^^^ Before using the EGO algorithm, we must create an initial kriging. In order to do this, we must create a design of experiment which fills the space. In this situation, the `LHSExperiment` is a good place to start (but other design of experiments may allow to better fill the space). We use a uniform distribution in order to create a LHS design with 50 points. .. GENERATED FROM PYTHON SOURCE LINES 79-86 .. code-block:: default listUniformDistributions = [ot.Uniform(lowerbound[i], upperbound[i]) for i in range(dim)] distribution = ot.ComposedDistribution(listUniformDistributions) sampleSize = 50 experiment = ot.LHSExperiment(distribution, sampleSize) inputSample = experiment.generate() outputSample = model(inputSample) .. GENERATED FROM PYTHON SOURCE LINES 87-92 .. code-block:: default graph = ot.Graph("Initial LHS design of experiment - n=%d" % (sampleSize), "$x_0$", "$x_1$", True) cloud = ot.Cloud(inputSample) graph.add(cloud) view = viewer.View(graph) .. image:: /auto_numerical_methods/optimization/images/sphx_glr_plot_ego_002.png :alt: Initial LHS design of experiment - n=50 :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 93-94 We now create the kriging metamodel. We selected the `SquaredExponential` covariance model with a constant basis (the `MaternModel` may perform better in this case). We use default settings (1.0) for the scale parameters of the covariance model, but configure the amplitude to 0.1, which better corresponds to the properties of the Ackley function. .. GENERATED FROM PYTHON SOURCE LINES 96-101 .. code-block:: default covarianceModel = ot.SquaredExponential([1.0] * dim, [0.5]) basis = ot.ConstantBasisFactory(dim).build() kriging = ot.KrigingAlgorithm(inputSample, outputSample, covarianceModel, basis) kriging.run() .. GENERATED FROM PYTHON SOURCE LINES 102-106 Create the optimization problem ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ We finally create the `OptimizationProblem` and solve it with `EfficientGlobalOptimization`. .. GENERATED FROM PYTHON SOURCE LINES 108-113 .. code-block:: default problem = ot.OptimizationProblem() problem.setObjective(model) bounds = ot.Interval(lowerbound, upperbound) problem.setBounds(bounds) .. GENERATED FROM PYTHON SOURCE LINES 114-117 In order to show the various options, we configure them all, even if most could be left to default settings in this case. The most important method is `setMaximumEvaluationNumber` which limits the number of iterations that the algorithm can perform. In the Ackley example, we choose to perform 10 iterations of the algorithm. .. GENERATED FROM PYTHON SOURCE LINES 119-124 .. code-block:: default algo = ot.EfficientGlobalOptimization(problem, kriging.getResult()) algo.setMaximumEvaluationNumber(10) algo.run() result = algo.getResult() .. GENERATED FROM PYTHON SOURCE LINES 125-127 .. code-block:: default result.getIterationNumber() .. rst-class:: sphx-glr-script-out Out: .. code-block:: none 10 .. GENERATED FROM PYTHON SOURCE LINES 128-130 .. code-block:: default result.getOptimalPoint() .. raw:: html

[-0.0975796,0.839969]



.. GENERATED FROM PYTHON SOURCE LINES 131-133 .. code-block:: default result.getOptimalValue() .. raw:: html

[3.00508]



.. GENERATED FROM PYTHON SOURCE LINES 134-136 .. code-block:: default fexact .. raw:: html

[4.44089e-16]



.. GENERATED FROM PYTHON SOURCE LINES 137-138 Compared to the minimum function value, we see that the EGO algorithm provides solution which is not very accurate. However, the optimal point is in the neighbourhood of the exact solution, and this is quite an impressive success given the limited amount of function evaluations: only 60 evaluations for the initial DOE and 10 iterations of the EGO algorithm, for a total equal to 70 function evaluations. .. GENERATED FROM PYTHON SOURCE LINES 140-143 .. code-block:: default graph = result.drawOptimalValueHistory() view = viewer.View(graph) .. image:: /auto_numerical_methods/optimization/images/sphx_glr_plot_ego_003.png :alt: Optimal value history :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 144-146 .. code-block:: default inputHistory = result.getInputSample() .. GENERATED FROM PYTHON SOURCE LINES 147-160 .. code-block:: default graph = model.draw(lowerbound, upperbound, [100]*dim) graph.setLegends([""]) graph.setTitle("Ackley function. Initial : black bullet. Solution : green diamond.") cloud = ot.Cloud(inputSample) cloud.setPointStyle("bullet") cloud.setColor("black") graph.add(cloud) cloud = ot.Cloud(inputHistory) cloud.setPointStyle("diamond") cloud.setColor("forestgreen") graph.add(cloud) view = viewer.View(graph) .. image:: /auto_numerical_methods/optimization/images/sphx_glr_plot_ego_004.png :alt: Ackley function. Initial : black bullet. Solution : green diamond. :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 161-164 We see that the initial (black) points are dispersed in the whole domain, while the solution points are much closer to the solution. However, the final solution produced by the EGO algorithm is not very accurate. This is why we finalize the process by adding a local optimization step. .. GENERATED FROM PYTHON SOURCE LINES 166-171 .. code-block:: default algo2 = ot.NLopt(problem, 'LD_LBFGS') algo2.setStartingPoint(result.getOptimalPoint()) algo2.run() result = algo2.getResult() .. GENERATED FROM PYTHON SOURCE LINES 172-174 .. code-block:: default result.getOptimalPoint() .. raw:: html

[4.5981e-07,0.952166]



.. GENERATED FROM PYTHON SOURCE LINES 175-176 The corrected solution is now extremely accurate. .. GENERATED FROM PYTHON SOURCE LINES 178-181 .. code-block:: default graph = result.drawOptimalValueHistory() view = viewer.View(graph) .. image:: /auto_numerical_methods/optimization/images/sphx_glr_plot_ego_005.png :alt: Optimal value history :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 182-186 Branin test-case ---------------- We now take a look at the :ref:`Branin-Hoo` function. .. GENERATED FROM PYTHON SOURCE LINES 188-190 Define the problem ^^^^^^^^^^^^^^^^^^ .. GENERATED FROM PYTHON SOURCE LINES 192-193 The Branin model is defined in the usecases module in a data class `BraninModel` : .. GENERATED FROM PYTHON SOURCE LINES 193-196 .. code-block:: default from openturns.usecases import branin_function as branin_function bm = branin_function.BraninModel() .. GENERATED FROM PYTHON SOURCE LINES 197-198 We load the dimension, .. GENERATED FROM PYTHON SOURCE LINES 198-200 .. code-block:: default dim = bm.dim .. GENERATED FROM PYTHON SOURCE LINES 201-202 the domain boundaries, .. GENERATED FROM PYTHON SOURCE LINES 202-205 .. code-block:: default lowerbound = bm.lowerbound upperbound = bm.upperbound .. GENERATED FROM PYTHON SOURCE LINES 206-207 and we load the model function : .. GENERATED FROM PYTHON SOURCE LINES 207-210 .. code-block:: default model = bm.model objectiveFunction = model.getMarginal(0) .. GENERATED FROM PYTHON SOURCE LINES 211-212 We build a sample out of the three minima : .. GENERATED FROM PYTHON SOURCE LINES 212-214 .. code-block:: default xexact = ot.Sample([bm.xexact1, bm.xexact2, bm.xexact3]) .. GENERATED FROM PYTHON SOURCE LINES 215-216 The minimum value attained `fexact` is : .. GENERATED FROM PYTHON SOURCE LINES 216-219 .. code-block:: default fexact = objectiveFunction(xexact) fexact .. raw:: html
y0
0-1.04741
1-1.04741
2-1.04741


.. GENERATED FROM PYTHON SOURCE LINES 220-224 .. code-block:: default graph = objectiveFunction.draw(lowerbound, upperbound, [100]*dim) graph.setTitle("Branin function") view = viewer.View(graph) .. image:: /auto_numerical_methods/optimization/images/sphx_glr_plot_ego_006.png :alt: Branin function :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 225-226 The Branin function has three local minimas. .. GENERATED FROM PYTHON SOURCE LINES 228-230 Create the initial kriging ^^^^^^^^^^^^^^^^^^^^^^^^^^ .. GENERATED FROM PYTHON SOURCE LINES 232-239 .. code-block:: default distribution = ot.ComposedDistribution([ot.Uniform(0.0, 1.0)] * dim) sampleSize = 50 experiment = ot.LHSExperiment(distribution, sampleSize) inputSample = experiment.generate() modelEval = model(inputSample) outputSample = modelEval.getMarginal(0) .. GENERATED FROM PYTHON SOURCE LINES 240-245 .. code-block:: default graph = ot.Graph("Initial LHS design of experiment - n=%d" % (sampleSize), "$x_0$", "$x_1$", True) cloud = ot.Cloud(inputSample) graph.add(cloud) view = viewer.View(graph) .. image:: /auto_numerical_methods/optimization/images/sphx_glr_plot_ego_007.png :alt: Initial LHS design of experiment - n=50 :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 246-250 .. code-block:: default covarianceModel = ot.SquaredExponential([1.0] * dim, [1.0]) basis = ot.ConstantBasisFactory(dim).build() kriging = ot.KrigingAlgorithm(inputSample, outputSample, covarianceModel, basis) .. GENERATED FROM PYTHON SOURCE LINES 251-255 .. code-block:: default noise = [x[1] for x in modelEval] kriging.setNoise(noise) kriging.run() .. GENERATED FROM PYTHON SOURCE LINES 256-258 Create and solve the problem ^^^^^^^^^^^^^^^^^^^^^^^^^^^^ .. GENERATED FROM PYTHON SOURCE LINES 260-261 We define the problem : .. GENERATED FROM PYTHON SOURCE LINES 261-266 .. code-block:: default problem = ot.OptimizationProblem() problem.setObjective(model) bounds = ot.Interval(lowerbound, upperbound) problem.setBounds(bounds) .. GENERATED FROM PYTHON SOURCE LINES 267-268 We configure the maximum number of function evaluations to 20. We assume that the function is noisy, with a constant variance. .. GENERATED FROM PYTHON SOURCE LINES 270-271 We configure the algorithm : .. GENERATED FROM PYTHON SOURCE LINES 271-280 .. code-block:: default algo = ot.EfficientGlobalOptimization(problem, kriging.getResult()) # assume constant noise var guessedNoiseFunction = 0.1 noiseModel = ot.SymbolicFunction(['x1', 'x2'], [str(guessedNoiseFunction)]) algo.setNoiseModel(noiseModel) algo.setMaximumEvaluationNumber(20) algo.run() result = algo.getResult() .. GENERATED FROM PYTHON SOURCE LINES 281-283 .. code-block:: default result.getIterationNumber() .. rst-class:: sphx-glr-script-out Out: .. code-block:: none 20 .. GENERATED FROM PYTHON SOURCE LINES 284-286 .. code-block:: default result.getOptimalPoint() .. raw:: html

[0.124428,0.803646]



.. GENERATED FROM PYTHON SOURCE LINES 287-289 .. code-block:: default result.getOptimalValue() .. raw:: html

[-1.04663]



.. GENERATED FROM PYTHON SOURCE LINES 290-292 .. code-block:: default fexact .. raw:: html
y0
0-1.04741
1-1.04741
2-1.04741


.. GENERATED FROM PYTHON SOURCE LINES 293-295 .. code-block:: default inputHistory = result.getInputSample() .. GENERATED FROM PYTHON SOURCE LINES 296-309 .. code-block:: default graph = objectiveFunction.draw(lowerbound, upperbound, [100]*dim) graph.setLegends([""]) graph.setTitle("Branin function. Initial : black bullet. Solution : green diamond.") cloud = ot.Cloud(inputSample) cloud.setPointStyle("bullet") cloud.setColor("black") graph.add(cloud) cloud = ot.Cloud(inputHistory) cloud.setPointStyle("diamond") cloud.setColor("forestgreen") graph.add(cloud) view = viewer.View(graph) .. image:: /auto_numerical_methods/optimization/images/sphx_glr_plot_ego_008.png :alt: Branin function. Initial : black bullet. Solution : green diamond. :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 310-311 We see that the EGO algorithm found the second local minimum. Given the limited number of function evaluations, the other local minimas have not been explored. .. GENERATED FROM PYTHON SOURCE LINES 313-317 .. code-block:: default graph = result.drawOptimalValueHistory() view = viewer.View(graph, axes_kw={"xticks": range(0, result.getIterationNumber(), 5)}) plt.show() .. image:: /auto_numerical_methods/optimization/images/sphx_glr_plot_ego_009.png :alt: Optimal value history :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 0 minutes 32.461 seconds) .. _sphx_glr_download_auto_numerical_methods_optimization_plot_ego.py: .. only :: html .. container:: sphx-glr-footer :class: sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_ego.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_ego.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_