6

scipy.optimize presents many different methods for local and global optimization of multivariate systems. However, I have a very long optimization run needed that may be interrupted (and in some cases I may want to interrupt it deliberately). Is there any way to restart... well, any of them? I mean, clearly one can provide the last, most optimized set of parameters found as the initial guess, but that's not the only parameter in play - for example, there are also gradients (jacobians, for example), populations in differential evolution, etc. I obviously don't want these to have to start over as well.

I see little way to prove these to scipy, nor to save its state. For functions that take a jacobian for example, there is a jacobian argument ("jac"), but it's either a boolean (indicating that your evaluation function returns a jacobian, which mine doesn't), or a callable function (I would only have the single result from the last run to provide). Nothing takes just an array of the last jacobian available. And with differential evolution, loss of the population would be horrible for performance and convergence.

Are there any solutions to this? Any way to resume optimizations at all?

KarenRei
  • 589
  • 6
  • 13
  • 1
    The best I've come up with is almost unthinkably ugly: fork the process to create one or more "breakpoint" versions (a process with an exact copy of your currently running optimization function - and whole application at that!) and have the breakpoint-versions hang in the optimizers' callback or evaluation calls until signaled to continue onward by some sort of IPC, or killed if rendered obsolete. Of course you couldn't change parameters between "runs" because they're still the same run! Also wouldn't persist through reboot. And just so ugly... – KarenRei Oct 12 '16 at 15:21
  • 1) if `jac` is your function, you can easily tell it to load and return `xx.jac` on the first call, with a *kwarg. 2) which optimizer specifically, differentialevolution ? 3) Do you want a tree of runs, which you want to run for a while, print/plot/save individually ? A large order; a sample sequence of commands: "A10 = run 10; A100 = run to 100; B = (A; params = ...; run 100) ..." woud clarify. – denis Oct 14 '16 at 11:43
  • I haven't yet benchmarked the different optimization methods to be able to know which one will work out best for me; I've still been tweaking my evaluation function for performance and accuracy. I have no clue how much for example it will tend to get trapped in local minima. "jac" is not my evaluation function, it's an argument that minimization functions take to specify how the jacobian (first order derivatives/gradient) is determined. A tree of runs would be nice but even being able to maintain just a single checkpoint would be great. How can you keep one and "run 100" iterations in scipy? – KarenRei Oct 14 '16 at 16:32
  • I've been trying to design my parameters so that crossover can be conducted in a way that doesn't produce nonsensical results (aka, trying to create parameters that are as independent from each other as possible) in case I need to use differential evolution. But differential evolution tends to be slower than gradient methods where local minima aren't a problem. So again it depends on how well the gradient methods converge. I'd also love to be able to run in parallel, but scipy doesn't seem to support that either... – KarenRei Oct 14 '16 at 16:38
  • 1
    For the record (just ran into this old thread), I found a painfully easy way to resume: cache all of the results of the minimization function, with the function parameters as their key. Save the cache and load it on restart. Add a check to the start of the minimization function to see if it's been run before. And that's it, done. Given the same parameters and the same returned results, minimization and differential evolution will reenact the exact same steps that they did before, only at lightning speed due to the cache. :) – KarenRei Jul 10 '17 at 20:35

2 Answers2

4

The general answer is no, there's no general solution apart from, just as you say, starting from the last estimate from the previous run.

For differential evolution specifically though, you can can instantiate the DifferentialEvolutionSolver, which you can pickle at checkpoint and unpickle to resume. (The suggestion comes from https://github.com/scipy/scipy/issues/6517)

ev-br
  • 24,968
  • 9
  • 65
  • 78
0

The following can save and restart from a previous x, but I gather you want to save and restart more state, e.g. gradients, too; can you clarify that ?

See also basinhopping, which has a nice-looking gui, pele-python .

#!/usr/bin/env python
""" Funcgradmn: wrap f() and grad(), save all x[] f[] grad[] to plot or restart """

from __future__ import division
import numpy as np

__version__ = "2016-10-18 oct denis"


class Funcgradmon(object):
    """ Funcgradmn: wrap f() and grad(), save all x[] f[] grad[] to plot or restart

    Example: minimize, save, restart --

    fg = Funcgradmon( func, gradfunc, verbose=1 )
        # fg(x): f(x), g(x)  for minimize( jac=True )

        # run 100 iter (if linesearch, 200-300 calls of fg()) --
    options = dict( maxiter=100 )  # ... 
    min0 = minimize( fg, x0, jac=True, options=options )
    fg.savez( "0.npz", paramstr="..." )  # to plot or restart

        # restart from x[50] --
        # (won't repeat the previous path from 50
        # unless you save and restore the whole state of the optimizer)
    x0 = fg.restart( 50 )
    # change params ...
    min50 = minimize( fg, x0, jac=True, options=options )
    """

    def __init__( self, func, gradfunc, verbose=1 ):
        self.func = func
        self.gradfunc = gradfunc
        self.verbose = verbose
        self.x, self.f, self.g = [], [], []  # growing lists
        self.t = 0

    def __call__( self, x ):
        """ f, g = func(x), gradfunc(x); save them; return f, g """
        x = np.asarray_chkfinite( x )  # always
        f = self.func(x)
        g = self.gradfunc(x)
        g = np.asarray_chkfinite( g )
        self.x.append( np.copy(x) )
        self.f.append( _copy( f ))
        self.g.append( np.copy(g) )
        if self.verbose:
            print "%3d:" % self.t ,
            fmt = "%-12g" if np.isscalar(f)  else "%s\t"
            print fmt % f ,
            print "x: %s" % x ,  # with user's np.set_printoptions
            print "\tgrad: %s" % g
                # better df dx dg
        # callback: plot
        self.t += 1
        return f, g

    def restart( self, n ):
        """ x0 = fg.restart( n )  returns x[n] to minimize( fg, x0 )
        """
        x0 = self.x[n]  # minimize from here
        del self.x[:n]
        del self.f[:n]
        del self.g[:n]
        self.t = n
        if self.verbose:
            print "Funcgradmon: restart from x[%d] %s" % (n, x0)
        return x0

    def savez( self, npzfile, **kw ):
        """ np.savez( npzfile, x= f= g= ) """
        x, f, g = map( np.array, [self.x, self.f, self.g] )
        if self.verbose:
            asum = "f: %s \nx: %s \ng: %s" % (
                _asum(f), _asum(x), _asum(g) )
            print "Funcgradmon: saving to %s: \n%s \n" % (npzfile, asum)
        np.savez( npzfile, x=x, f=f, g=g, **kw )

    def load( self, npzfile ):
        load = np.load( npzfile )
        x, f, g = load["x"], load["f"], load["g"]
        if self.verbose:
            asum = "f: %s \nx: %s \ng: %s" % (
                _asum(f), _asum(x), _asum(g) )
            print "Funcgradmon: load %s: \n%s \n" % (npzfile, asum)
        self.x = list( x )
        self.f = list( f )
        self.g = list( g )
        self.loaddict = load
        return self.restart( len(x) - 1 )


def _asum( X ):
    """ one-line array summary: "shape type min av max" """
    if not hasattr( X, "dtype" ):
        return str(X)
    return "%s %s  min av max %.3g %.3g %.3g" % (
            X.shape, X.dtype, X.min(), X.mean(), X.max() )

def _copy( x ):
    return x if x is None  or np.isscalar(x) \
        else np.copy( x )

#...............................................................................
if __name__ == "__main__":
    import sys
    from scipy.optimize import minimize, rosen, rosen_der

    np.set_printoptions( threshold=20, edgeitems=10, linewidth=140,
            formatter = dict( float = lambda x: "%.3g" % x ))  # float arrays %.3g

    dim = 3
    method = "cg"
    maxiter = 10  # 1 linesearch -> 2-3 calls of fg

    # to change these params, run this.py a=1 b=None 'c = ...'  in sh or ipython
    for arg in sys.argv[1:]:
        exec( arg )

    print "\n", 80 * "-"
    print "Funcgradmon: dim %d  method %s  maxiter %d \n" % (
            dim, method, maxiter )
    x0 = np.zeros( dim )

    #...........................................................................
    fg = Funcgradmon( rosen, rosen_der, verbose=1 )
    options = dict( maxiter=maxiter )  # ... 

    min0 = minimize( fg, x0, jac=True, method=method, options=options )
    fg.savez( "0.npz", paramstr="..." )  # to plot or restart

    x0 = fg.restart( 5 )  # = fg.x[5]
    # change params, print them all
    min5 = minimize( fg, x0, jac=True, method=method, options=options )

    fg.savez( "5.npz", paramstr="..." )
denis
  • 21,378
  • 10
  • 65
  • 88