Squashed commit of #54 - Add Optimisers

examples/scripts/rmse_estimation.py

@@ -45,7 +45,7 @@
 )
 
 # Run the optimisation problem
-results, last_optim, num_evals = parameterisation.run()
+x, output, final_cost, num_evals = parameterisation.run()
 
 
 # get MAP estimate, starting at a random initial point in parameter space

pybop/__init__.py

@@ -1,13 +1,13 @@
 #
-# Root of the pybop module.
+# Root of the PyBOP module.
 # Provides access to all shared functionality (models, solvers, etc.).
 #
 # This file is adapted from Pints
 # (see https://github.com/pints-team/pints)
 #
 
 import sys
-import os
+from os import path
 
 #
 # Version info
@@ -20,8 +20,8 @@
 # Float format: a float can be converted to a 17 digit decimal and back without
 # loss of information
 FLOAT_FORMAT = "{: .17e}"
-# Absolute path to the pybop module
-script_path = os.path.dirname(__file__)
+# Absolute path to the PyBOP repo
+script_path = path.dirname(__file__)
 
 #
 # Cost function class
@@ -50,6 +50,7 @@
 from .optimisers.base_optimiser import BaseOptimiser
 from .optimisers.nlopt_optimize import NLoptOptimize
 from .optimisers.scipy_minimize import SciPyMinimize
+from .optimisers.pints_optimiser import PintsOptimiser, PintsError, PintsBoundaries
 
 #
 # Parameter classes
@@ -64,6 +65,6 @@
 from .plotting.quick_plot import QuickPlot
 
 #
-# Remove any imported modules, so we don't expose them as part of pybop
+# Remove any imported modules, so we don't expose them as part of PyBOP
 #
 del sys

pybop/optimisers/base_optimiser.py

@@ -8,14 +8,13 @@ class BaseOptimiser:
     def __init__(self):
         self.name = "Base Optimiser"
 
-    def optimise(self, cost_function, x0, bounds, method=None):
+    def optimise(self, cost_function, x0, bounds):
         """
         Optimisiation method to be overloaded by child classes.
 
         """
         self.cost_function = cost_function
         self.x0 = x0
-        self.method = method
         self.bounds = bounds
 
         # Run optimisation

pybop/optimisers/nlopt_optimize.py

@@ -7,36 +7,43 @@ class NLoptOptimize(BaseOptimiser):
     Wrapper class for the NLOpt optimiser class. Extends the BaseOptimiser class.
     """
 
-    def __init__(self, method=None, x0=None, xtol=None):
+    def __init__(self, x0, xtol=None, method=None):
         super().__init__()
         self.name = "NLOpt Optimiser"
 
         if method is not None:
-            self.opt = nlopt.opt(method, len(x0))
+            self.optim = nlopt.opt(method, len(x0))
         else:
-            self.opt = nlopt.opt(nlopt.LN_BOBYQA, len(x0))
+            self.optim = nlopt.opt(nlopt.LN_BOBYQA, len(x0))
 
         if xtol is not None:
-            self.opt.set_xtol_rel(xtol)
+            self.optim.set_xtol_rel(xtol)
         else:
-            self.opt.set_xtol_rel(1e-5)
+            self.optim.set_xtol_rel(1e-5)
 
     def _runoptimise(self, cost_function, x0, bounds):
         """
-        Run the NLOpt opt method.
+        Run the NLOpt optimisation method.
 
-        Parameters
+        Inputs
         ----------
         cost_function: function for optimising
-        method: optimisation method
-        x0: Initialisation array
+        method: optimisation algorithm
+        x0: initialisation array
         bounds: bounds array
         """
 
-        self.opt.set_min_objective(cost_function)
-        self.opt.set_lower_bounds(bounds["lower"])
-        self.opt.set_upper_bounds(bounds["upper"])
-        results = self.opt.optimize(x0)
-        num_evals = self.opt.get_numevals()
+        # Pass settings to the optimiser
+        self.optim.set_min_objective(cost_function)
+        self.optim.set_lower_bounds(bounds["lower"])
+        self.optim.set_upper_bounds(bounds["upper"])
 
-        return results, self.opt.last_optimum_value(), num_evals
+        # Run the optimser
+        x = self.optim.optimize(x0)
+
+        # Get performance statistics
+        output = self.optim
+        final_cost = self.optim.last_optimum_value()
+        num_evals = self.optim.get_numevals()
+
+        return x, output, final_cost, num_evals

pybop/optimisers/pints_optimiser.py

@@ -0,0 +1,158 @@
+import pybop
+import pints
+from pybop.optimisers.base_optimiser import BaseOptimiser
+from pints import ErrorMeasure
+
+
+class PintsOptimiser(BaseOptimiser):
+    """
+    Wrapper class for the PINTS optimisation class. Extends the BaseOptimiser class.
+    """
+
+    def __init__(self, x0, xtol=None, method=None):
+        super().__init__()
+        self.name = "PINTS Optimiser"
+
+        if method is not None:
+            self.method = method
+        else:
+            self.method = pints.PSO
+
+    def _runoptimise(self, cost_function, x0, bounds):
+        """
+        Run the PINTS optimisation method.
+
+        Inputs
+        ----------
+        cost_function: function for optimising
+        method: optimisation algorithm
+        x0: initialisation array
+        bounds: bounds array
+        """
+
+        # Wrap bounds
+        boundaries = pybop.PintsBoundaries(bounds, x0)
+
+        # Wrap error measure
+        error = pybop.PintsError(cost_function, x0)
+
+        # Set up optimisation controller
+        controller = pints.OptimisationController(
+            error, x0, boundaries=boundaries, method=self.method
+        )
+        controller.set_max_unchanged_iterations(20)  # default 200
+
+        # Run the optimser
+        x, final_cost = controller.run()
+
+        # Get performance statistics
+        # output = *pass all output*
+        # final_cost
+        # num_evals
+        output = None
+        num_evals = None
+
+        return x, output, final_cost, num_evals
+
+
+class PintsError(ErrorMeasure):
+    """
+    An interface class for PyBOP that extends the PINTS ErrorMeasure class.
+
+    From PINTS:
+    Abstract base class for objects that calculate some scalar measure of
+    goodness-of-fit (for a model and a data set), such that a smaller value
+    means a better fit.
+
+    ErrorMeasures are callable objects: If ``e`` is an instance of an
+    :class:`ErrorMeasure` class you can calculate the error by calling ``e(p)``
+    where ``p`` is a point in parameter space.
+    """
+
+    def __init__(self, cost_function, x0):
+        self.cost_function = cost_function
+        self.x0 = x0
+
+    def __call__(self, x):
+        cost = self.cost_function(x)
+
+        return cost
+
+    def evaluateS1(self, x):
+        """
+        Evaluates this error measure, and returns the result plus the partial
+        derivatives of the result with respect to the parameters.
+
+        The returned data has the shape ``(e, e')`` where ``e`` is a scalar
+        value and ``e'`` is a sequence of length ``n_parameters``.
+
+        *This is an optional method that is not always implemented.*
+        """
+        raise NotImplementedError
+
+    def n_parameters(self):
+        """
+        Returns the dimension of the parameter space this measure is defined
+        over.
+        """
+        return len(self.x0)
+
+
+class PintsBoundaries(object):
+    """
+    An interface class for PyBOP that extends the PINTS ErrorMeasure class.
+
+    From PINTS:
+    Abstract class representing boundaries on a parameter space.
+    """
+
+    def __init__(self, bounds, x0):
+        self.bounds = bounds
+        self.x0 = x0
+
+    def check(self, parameters):
+        """
+        Returns ``True`` if and only if the given point in parameter space is
+        within the boundaries.
+
+        Parameters
+        ----------
+        parameters
+            A point in parameter space
+        """
+        result = False
+        if (
+            parameters[0] >= self.bounds["lower"][0]
+            and parameters[1] >= self.bounds["lower"][1]
+            and parameters[0] <= self.bounds["upper"][0]
+            and parameters[1] <= self.bounds["upper"][1]
+        ):
+            result = True
+
+        return result
+
+    def n_parameters(self):
+        """
+        Returns the dimension of the parameter space these boundaries are
+        defined on.
+        """
+        return len(self.x0)
+
+    def sample(self, n=1):
+        """
+        Returns ``n`` random samples from within the boundaries, for example to
+        use as starting points for an optimisation.
+
+        The returned value is a NumPy array with shape ``(n, d)`` where ``n``
+        is the requested number of samples, and ``d`` is the dimension of the
+        parameter space these boundaries are defined on.
+
+        *Note that implementing :meth:`sample()` is optional, so some boundary
+        types may not support it.*
+
+        Parameters
+        ----------
+        n : int
+            The number of points to sample
+        """
+        raise NotImplementedError

pybop/optimisers/scipy_minimize.py

@@ -4,38 +4,50 @@
 
 class SciPyMinimize(BaseOptimiser):
     """
-    Wrapper class for the Scipy optimiser class. Extends the BaseOptimiser class.
+    Wrapper class for the SciPy optimisation class. Extends the BaseOptimiser class.
     """
 
-    def __init__(self, cost_function, x0, bounds=None, options=None):
+    def __init__(self, x0, xtol=None, method=None, options=None):
         super().__init__()
-        self.cost_function = cost_function
-        self.method = options.optmethod
-        self.x0 = x0 or cost_function.x0
-        self.bounds = bounds
-        self.options = options
         self.name = "Scipy Optimiser"
 
-    def _runoptimise(self):
+        if method is None:
+            self.method = method
+        else:
+            self.method = "BFGS"
+
+        if xtol is not None:
+            self.xtol = xtol
+        else:
+            self.xtol = 1e-5
+
+        self.options = options
+
+    def _runoptimise(self, cost_function, x0, bounds):
         """
-        Run the Scipy opt method.
+        Run the SciPy optimisation method.
 
-        Parameters
+        Inputs
         ----------
         cost_function: function for optimising
-        method: optimisation method
-        x0: Initialisation array
-        options: options dictionary
+        method: optimisation algorithm
+        x0: initialisation array
         bounds: bounds array
         """
 
-        if self.method is not None and self.bounds is not None:
-            opt = minimize(
-                self.cost_function, self.x0, method=self.method, bounds=self.bounds
-            )
-        elif self.method is not None:
-            opt = minimize(self.cost_function, self.x0, method=self.method)
-        else:
-            opt = minimize(self.cost_function, self.x0, method="BFGS")
+        # Reformat bounds
+        bounds = (
+            (lower, upper) for lower, upper in zip(bounds["lower"], bounds["upper"])
+        )
+
+        # Run the optimser
+        output = minimize(
+            cost_function, x0, method=self.method, bounds=bounds, tol=self.xtol
+        )
+
+        # Get performance statistics
+        x = output.x
+        final_cost = output.fun
+        num_evals = output.nfev
 
-        return opt
+        return x, output, final_cost, num_evals

pybop/parameters/base_parameter.py

@@ -1,6 +1,6 @@
 class Parameter:
     """ ""
-    Class for creating parameters in pybop.
+    Class for creating parameters in PyBOP.
     """
 
     def __init__(self, name, value=None, prior=None, bounds=None):

pybop/parameters/base_parameter_set.py

@@ -3,7 +3,7 @@
 
 class ParameterSet:
     """
-    Class for creating parameter sets in pybop.
+    Class for creating parameter sets in PyBOP.
     """
 
     def __new__(cls, method, name):