From 6594aab9bb1612d698add28dc8736cf788a39aa2 Mon Sep 17 00:00:00 2001 From: Mike McKerns Date: Thu, 27 Jun 2013 23:25:23 -0700 Subject: [PATCH 01/10] Initial commit --- README.md | 4 ++++ 1 file changed, 4 insertions(+) create mode 100644 README.md diff --git a/README.md b/README.md new file mode 100644 index 00000000..eac9bc26 --- /dev/null +++ b/README.md @@ -0,0 +1,4 @@ +mystic +====== + +a framework for highly-constrained non-convex optimization and uncertainty quantification From 77aed4090f3c41ea935dd1d9007a4e5b3cfeabbb Mon Sep 17 00:00:00 2001 From: uber Date: Fri, 28 Jun 2013 16:05:47 -0500 Subject: [PATCH 02/10] error found with no scipy read_array method --- examples/test_br8.py | 13 +++++++++---- 1 file changed, 9 insertions(+), 4 deletions(-) diff --git a/examples/test_br8.py b/examples/test_br8.py index 34f21b84..fb771319 100755 --- a/examples/test_br8.py +++ b/examples/test_br8.py @@ -121,10 +121,15 @@ def curry(x): # 2D (a fine mesh solution can be computed by test_br8_mpi.py) try: - import scipy.io, pylab - X = scipy.io.read_array(open('test_br8_mpi.out.X')) - Y = scipy.io.read_array(open('test_br8_mpi.out.Y')) - V = scipy.io.read_array(open('test_br8_mpi.out.V')) + # found error scipy no longer uses this method + #import scipy.io, pylab + import numpy, pylab + #X = scipy.io.read_array(open('test_br8_mpi.out.X')) + #Y = scipy.io.read_array(open('test_br8_mpi.out.Y')) + #V = scipy.io.read_array(open('test_br8_mpi.out.V')) + X = numpy.loadtxt('test_br8_mpi.out.X') + Y = numpy.loadtxt('test_br8_mpi.out.Y') + V = numpy.loadtxt('test_br8_mpi.out.V') pylab.clf() pylab.plot([[a4]],[[a5]],'k+') pylab.xlabel('a4') From 11a0ea1d3ed0f2410e8f57a945236cd90f8ee460 Mon Sep 17 00:00:00 2001 From: uber Date: Fri, 28 Jun 2013 16:25:20 -0500 Subject: [PATCH 03/10] issue-1: bug fixed --- examples/test_br8.py | 5 ----- 1 file changed, 5 deletions(-) diff --git a/examples/test_br8.py b/examples/test_br8.py index fb771319..a18a4f3b 100755 --- a/examples/test_br8.py +++ b/examples/test_br8.py @@ -121,12 +121,7 @@ def curry(x): # 2D (a fine mesh solution can be computed by test_br8_mpi.py) try: - # found error scipy no longer uses this method - #import scipy.io, pylab import numpy, pylab - #X = scipy.io.read_array(open('test_br8_mpi.out.X')) - #Y = scipy.io.read_array(open('test_br8_mpi.out.Y')) - #V = scipy.io.read_array(open('test_br8_mpi.out.V')) X = numpy.loadtxt('test_br8_mpi.out.X') Y = numpy.loadtxt('test_br8_mpi.out.Y') V = numpy.loadtxt('test_br8_mpi.out.V') From 45186f1df294f4ea3e057ea8c912a684d89a8a09 Mon Sep 17 00:00:00 2001 From: Mike McKerns Date: Thu, 27 Jun 2013 23:25:23 -0700 Subject: [PATCH 04/10] Initial commit --- README.md | 4 ++++ 1 file changed, 4 insertions(+) create mode 100644 README.md diff --git a/README.md b/README.md new file mode 100644 index 00000000..eac9bc26 --- /dev/null +++ b/README.md @@ -0,0 +1,4 @@ +mystic +====== + +a framework for highly-constrained non-convex optimization and uncertainty quantification From 87ddd469ddd690a7db84de5922b2640d67b054c0 Mon Sep 17 00:00:00 2001 From: Mike McKerns Date: Thu, 11 Jul 2013 10:21:14 -0700 Subject: [PATCH 05/10] merged in updates to README.md from svn --- README.md | 88 ++++++++++++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 87 insertions(+), 1 deletion(-) diff --git a/README.md b/README.md index eac9bc26..36d326a1 100644 --- a/README.md +++ b/README.md @@ -1,4 +1,90 @@ mystic ====== - a framework for highly-constrained non-convex optimization and uncertainty quantification + +About Mystic +------------ +The mystic framework provides a collection of optimization algorithms +and tools that allows the user to more robustly (and readily) solve +optimization problems. All optimization algorithms included in mystic +provide workflow at the fitting layer, not just access to the algorithms +as function calls. Mystic gives the user fine-grained power to both +monitor and steer optimizations as the fit processes are running. + +Where possible, mystic optimizers share a common interface, and thus can +be easily swapped without the user having to write any new code. Mystic +solvers all conform to a solver API, thus also have common method calls +to configure and launch an optimization job. For more details, see +`mystic.abstract_solver`. The API also makes it easy to bind a favorite +3rd party solver into the mystic framework. + +By providing a robust interface designed to allow the user to easily +configure and control solvers, mystic reduces the barrier to implementing +a target fitting problem as stable code. Thus the user can focus on +building their physical models, and not spend time hacking together an +interface to optimization code. + +Mystic is in the early development stages, and any user feedback is +highly appreciated. Contact Mike McKerns [mmckerns at caltech dot edu] +with comments, suggestions, and any bugs you may find. A list of known +issues is maintained at http://trac.mystic.cacr.caltech.edu/project/mystic/query. + +Major Features +-------------- +Mystic provides a stock set of configurable, controllable solvers with:: + * a common interface + * the ability to impose solver-independent bounds constraints + * the ability to apply solver-independent monitors + * the ability to configure solver-independent termination conditions + * a control handler yielding: [pause, continue, exit, and user_callback] + * ease in selecting initial conditions: [initial_guess, random] + * ease in selecting mutation strategies (for differential evolution) + +To get up and running quickly, mystic also provides infrastructure to:: + * easily generate a fit model (several example models are included) + * configure and auto-generate a cost function from a model + * extend fit jobs to parallel & distributed resources + * couple models with optimization parameter constraints [COMING SOON] + + +Current Release +--------------- +The latest released version of mystic is available from:: + http://trac.mystic.cacr.caltech.edu/project/mystic + +Mystic is distributed under a modified BSD license. + +Development Release +------------------- +You can get the latest development release with all the shiny new features at:: + http://dev.danse.us/packages. + +or even better, fork us on our github mirror of the svn trunk:: + https://github.com/uqfoundation + +Citation +-------- +If you use mystic to do research that leads to publication, we ask that you +acknowledge use of mystic by citing the following in your publication:: + + M.M. McKerns, L. Strand, T. Sullivan, A. Fang, M.A.G. Aivazis, + "Building a framework for predictive science", Proceedings of + the 10th Python in Science Conference, 2011; + http://arxiv.org/pdf/1202.1056 + + Michael McKerns, Patrick Hung, and Michael Aivazis, + "mystic: a simple model-independent inversion framework", 2009- ; + http://trac.mystic.cacr.caltech.edu/project/mystic + +More Information +---------------- +Probably the best way to get started is to look at the tutorial examples provided +within the user's guide. The source code is also generally well documented, +so further questions may be resolved by inspecting the code itself, or through +browsing the reference manual. For those who like to leap before +they look, you can jump right to the installation instructions. If the aforementioned documents +do not adequately address your needs, please send us feedback. + +Mystic is an active research tool. There are a growing number of publications and presentations that +discuss real-world examples and new features of mystic in greater detail than presented in the user's guide. +If you would like to share how you use mystic in your work, please send us a link. From d12f3a39abf433c305b9e368f5d555227d7d9584 Mon Sep 17 00:00:00 2001 From: Mike McKerns Date: Thu, 11 Jul 2013 11:20:26 -0700 Subject: [PATCH 06/10] fixed formatting in README.md --- README.md | 24 +++++++++++++----------- 1 file changed, 13 insertions(+), 11 deletions(-) diff --git a/README.md b/README.md index 36d326a1..24748f2a 100644 --- a/README.md +++ b/README.md @@ -32,19 +32,21 @@ issues is maintained at http://trac.mystic.cacr.caltech.edu/project/mystic/query Major Features -------------- Mystic provides a stock set of configurable, controllable solvers with:: - * a common interface - * the ability to impose solver-independent bounds constraints - * the ability to apply solver-independent monitors - * the ability to configure solver-independent termination conditions - * a control handler yielding: [pause, continue, exit, and user_callback] - * ease in selecting initial conditions: [initial_guess, random] - * ease in selecting mutation strategies (for differential evolution) + +* a common interface +* the ability to impose solver-independent bounds constraints +* the ability to apply solver-independent monitors +* the ability to configure solver-independent termination conditions +* a control handler yielding: [pause, continue, exit, and user_callback] +* ease in selecting initial conditions: [initial_guess, random] +* ease in selecting mutation strategies (for differential evolution) To get up and running quickly, mystic also provides infrastructure to:: - * easily generate a fit model (several example models are included) - * configure and auto-generate a cost function from a model - * extend fit jobs to parallel & distributed resources - * couple models with optimization parameter constraints [COMING SOON] + +* easily generate a fit model (several example models are included) +* configure and auto-generate a cost function from a model +* extend fit jobs to parallel & distributed resources +* couple models with optimization parameter constraints [COMING SOON] Current Release From bf23b0e8188439b1cb09a776519a64076dbd83e5 Mon Sep 17 00:00:00 2001 From: Erik Welch Date: Sat, 31 May 2014 12:52:33 -0400 Subject: [PATCH 07/10] Fix typo from 8995009: missing apostrophe in setup.py --- setup.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/setup.py b/setup.py index 1f56dd84..8593de55 100644 --- a/setup.py +++ b/setup.py @@ -304,7 +304,7 @@ def write_info_py(filename='mystic/info.py'): numpy_version = '>=1.0' sympy_version = '>=0.6.7' dill_version = '>=0.2' -klepto_version = '>=0.1 +klepto_version = '>=0.1' scipy_version = '>=0.6.0' matplotlib_version = '>=0.91' if has_setuptools: From d29acb56ea6492c5592e773d00903abe0658e7f5 Mon Sep 17 00:00:00 2001 From: Christoph Deil Date: Sat, 4 Oct 2014 23:10:42 +0200 Subject: [PATCH 08/10] Add .gitignore file --- .gitignore | 11 +++++++++++ 1 file changed, 11 insertions(+) create mode 100644 .gitignore diff --git a/.gitignore b/.gitignore new file mode 100644 index 00000000..b05b2901 --- /dev/null +++ b/.gitignore @@ -0,0 +1,11 @@ +*.pyc +__pycache__ +.project +.pydevproject +.settings +.idea +build +dist +mystic.egg-info +mystic/info.py +README From 223cba11dedee0ec3657201cabd1d866469b5176 Mon Sep 17 00:00:00 2001 From: Christoph Deil Date: Sat, 4 Oct 2014 23:17:00 +0200 Subject: [PATCH 09/10] Change comparison to None to get rid of FutureWarning --- mystic/abstract_solver.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/mystic/abstract_solver.py b/mystic/abstract_solver.py index 4594cf0a..dc084873 100755 --- a/mystic/abstract_solver.py +++ b/mystic/abstract_solver.py @@ -377,16 +377,16 @@ def SetRandomInitialPoints(self, min=None, max=None): input:: - min, max: must be a sequence of length self.nDim - each min[i] should be <= the corresponding max[i]""" - if min == None: min = self._defaultMin - if max == None: max = self._defaultMax + if min is None: min = self._defaultMin + if max is None: max = self._defaultMax #if numpy.any(( asarray(min) > asarray(max) ),0): # raise ValueError, "each min[i] must be <= the corresponding max[i]" if len(min) != self.nDim or len(max) != self.nDim: raise ValueError, "bounds array must be length %s" % self.nDim # when 'some' of the bounds are given as 'None', replace with default for i in range(len(min)): - if min[i] == None: min[i] = self._defaultMin[0] - if max[i] == None: max[i] = self._defaultMax[0] + if min[i] is None: min[i] = self._defaultMin[0] + if max[i] is None: max[i] = self._defaultMax[0] import random #generate random initial values for i in range(len(self.population)): From 9903413f4cad733070ef00c6337364b2947f00dd Mon Sep 17 00:00:00 2001 From: Jean-Christophe Fillion-Robin Date: Wed, 8 Jul 2015 03:16:02 -0400 Subject: [PATCH 10/10] Re-factor model_plotter for better integration with IPython This commit implements the following improvements: (1) adds "mystic.model_plotter.model_plotter(cmdargs)" function (2) adds an additional parameter -p/--plot-filepath 1) mystic.model_plotter.model_plotter(cmdargs) This convenience functions allows to directly invoke the plotter providing either (a) a list of arguments or (b) a string representing the arguments. For example, from IPython, ones could plot using the following approaches: %matplotlib inline from mystic.model_plotter import model_plotter model_plotter(cmdargs='mystic.models.zimmermann log.txt -b "-5:10:.1, -5:10:.1" -d -x 1') or %matplotlib inline from mystic.model_plotter import model_plotter model_plotter(cmdargs=['mystic.models.zimmermann', 'log.txt', '-b', '-5:10:.1, -5:10:.1', '-d', '-x', '1']) 2) additional parameter -p/--plot-filepath This new argument allows to specify an image filename where the generated plot should be saved. For example: mystic_model_plotter.py mystic.models.zimmermann -b "-5:10:.1, -5:10:.1" -d -x 1 -p output.png --- mystic/model_plotter.py | 584 ++++++++++++++++++++++++++++++++ scripts/mystic_model_plotter.py | 539 +---------------------------- 2 files changed, 589 insertions(+), 534 deletions(-) create mode 100644 mystic/model_plotter.py diff --git a/mystic/model_plotter.py b/mystic/model_plotter.py new file mode 100644 index 00000000..050d83d9 --- /dev/null +++ b/mystic/model_plotter.py @@ -0,0 +1,584 @@ +#!/usr/bin/env python +# +# Author: Patrick Hung (patrickh @caltech) +# Author: Jean-Christophe Fillion-Robin (jchris.fillionr @kitware.com) +# Copyright (c) 1997-2015 California Institute of Technology. +# License: 3-clause BSD. The full license text is available at: +# - http://trac.mystic.cacr.caltech.edu/project/mystic/browser/mystic/LICENSE + +__doc__ = """ +mystic_model_plotter.py [options] model (filename) + +generate surface contour plots for model, specified by full import path +generate model trajectory from logfile (or solver restart file), if provided + +The option "bounds" takes an indicator string, where the bounds should +be given as comma-separated slices. For example, using bounds = "-1:10, 0:20" +will set the lower and upper bounds for x to be (-1,10) and y to be (0,20). +The "step" can also be given, to control the number of lines plotted in the +grid. Thus "-1:10:.1, 0:20" would set the bounds as above, but use increments +of .1 along x and the default step along y. For models with > 2D, the bounds +can be used to specify 2 dimensions plus fixed values for remaining dimensions. +Thus, "-1:10, 0:20, 1.0" would plot the 2D surface where the z-axis was fixed +at z=1.0. + +The option "label" takes comma-separated strings. For example, label = "x,y," +will place 'x' on the x-axis, 'y' on the y-axis, and nothing on the z-axis. +LaTeX is also accepted. For example, label = "$ h $, $ {\\alpha}$, $ v$" will +label the axes with standard LaTeX math formatting. Note that the leading +space is required, while a trailing space aligns the text with the axis +instead of the plot frame. + +The option "reduce" can be given to reduce the output of a model to a scalar, +thus converting 'model(params)' to 'reduce(model(params))'. A reducer is given +by the import path (e.g. 'numpy.add'). The option "scale" will convert the plot +to log-scale, and scale the cost by 'z=log(4*z*scale+1)+2'. This is useful for +visualizing small contour changes around the minimium. If using log-scale +produces negative numbers, the option "shift" can be used to shift the cost +by 'z=z+shift'. Both shift and scale are intended to help visualize contours. + +Required Inputs: + model full import path for the model (e.g. mystic.models.rosen) + +Additional Inputs: + filename name of the convergence logfile (e.g. log.txt) +""" + +from mpl_toolkits.mplot3d import axes3d +import matplotlib.pyplot as plt +from matplotlib import cm + +from mystic.munge import read_history +from mystic.munge import logfile_reader, raw_to_support + +#XXX: better if reads single id only? (e.g. same interface as read_history) +def get_history(source, ids=None): + """get params and cost from the given source + +source is the name of the trajectory logfile (or solver instance) +if provided, ids are the list of 'run ids' to select + """ + try: # if it's a logfile, it might be multi-id + step, param, cost = logfile_reader(source) + except: # it's not a logfile, so read and return + param, cost = read_history(source) + return [param],[cost] + + # split (i,id) into iteration and id + multinode = len(step[0]) - 1 #XXX: what if step = []? + if multinode: id = [i[1] for i in step] + else: id = [0 for i in step] + + params = [[] for i in range(max(id) + 1)] + costs = [[] for i in range(len(params))] + # populate params for each id with the corresponding (param,cost) + for i in range(len(id)): + if ids is None or id[i] in ids: # take only the selected 'id' + params[id[i]].append(param[i]) + costs[id[i]].append(cost[i]) + params = [r for r in params if len(r)] # only keep selected 'ids' + costs = [r for r in costs if len(r)] # only keep selected 'ids' + + # convert to support format + for i in range(len(params)): + params[i], costs[i] = raw_to_support(params[i], costs[i]) + return params, costs + + +def get_instance(location, *args, **kwds): + """given the import location of a model or model class, return the model + +args and kwds will be passed to the constructor of the model class + """ + package, target = location.rsplit('.',1) + exec "from %s import %s as model" % (package, target) + import inspect + if inspect.isclass(model): + model = model(*args, **kwds) + return model + + +def parse_input(option): + """parse 'option' string into 'select', 'axes', and 'mask' + +select contains the dimension specifications on which to plot +axes holds the indicies of the parameters selected to plot +mask is a dictionary of the parameter indicies and fixed values + +For example, + >>> select, axes, mask = parse_input("-1:10:.1, 0.0, 5.0, -50:50:.5") + >>> select + [0, 3] + >>> axes + "-1:10:.1, -50:50:.5" + >>> mask + {1: 0.0, 2: 5.0} + """ + option = option.split(',') + select = [] + axes = [] + mask = {} + for index,value in enumerate(option): + if ":" in value: + select.append(index) + axes.append(value) + else: + mask.update({index:float(value)}) + axes = ','.join(axes) + return select, axes, mask + + +def parse_axes(option, grid=True): + """parse option string into grid axes; using modified numpy.ogrid notation + +For example: + option='-1:10:.1, 0:10:.1' yields x,y=ogrid[-1:10:.1,0:10:.1], + +If grid is False, accept options suitable for line plotting. +For example: + option='-1:10' yields x=ogrid[-1:10] and y=0, + option='-1:10, 2' yields x=ogrid[-1:10] and y=2, + +Returns tuple (x,y) with 'x,y' defined above. + """ + import numpy + option = option.split(',') + opt = dict(zip(['x','y','z'],option)) + if len(option) > 2 or len(option) < 1: + raise ValueError("invalid format string: '%s'" % ','.join(option)) + z = bool(grid) + if len(option) == 1: opt['y'] = '0' + xd = True if ':' in opt['x'] else False + yd = True if ':' in opt['y'] else False + #XXX: accepts option='3:1', '1:1', and '1:2:10' (try to catch?) + if xd and yd: + try: # x,y form a 2D grid + exec('x,y = numpy.ogrid[%s,%s]' % (opt['x'],opt['y'])) + except: # AttributeError: + raise ValueError("invalid format string: '%s'" % ','.join(option)) + elif xd and not z: + try: + exec('x = numpy.ogrid[%s]' % opt['x']) + y = float(opt['y']) + except: # (AttributeError, SyntaxError, ValueError): + raise ValueError("invalid format string: '%s'" % ','.join(option)) + elif yd and not z: + try: + x = float(opt['x']) + exec('y = numpy.ogrid[%s]' % opt['y']) + except: # (AttributeError, SyntaxError, ValueError): + raise ValueError("invalid format string: '%s'" % ','.join(option)) + else: + raise ValueError("invalid format string: '%s'" % ','.join(option)) + if not x.size or not y.size: + raise ValueError("invalid format string: '%s'" % ','.join(option)) + return x,y + + +def draw_projection(x, cost, scale=True, shift=False, style=None, figure=None): + """draw a solution trajectory (for overlay on a 1D plot) + +x is the sequence of values for one parameter (i.e. a parameter trajectory) +cost is the sequence of costs (i.e. the solution trajectory) +if scale is provided, scale the intensity as 'z = log(4*z*scale+1)+2' +if shift is provided, shift the intensity as 'z = z+shift' (useful for -z's) +if style is provided, set the line style (e.g. 'w-o', 'k-', 'ro') +if figure is provided, plot to an existing figure + """ + if not figure: figure = plt.figure() + ax = figure.gca() + ax.autoscale(tight=True) + + if style in [None, False]: + style = 'k-o' + import numpy + if shift: + if shift is True: #NOTE: MAY NOT be the exact minimum + shift = max(-numpy.min(cost), 0.0) + 0.5 # a good guess + cost = numpy.asarray(cost)+shift + cost = numpy.asarray(cost) + if scale: + cost = numpy.log(4*cost*scale+1)+2 + + ax.plot(x,cost, style, linewidth=2, markersize=4) + #XXX: need to 'correct' the z-axis (or provide easy conversion) + return figure + + +def draw_trajectory(x, y, cost=None, scale=True, shift=False, style=None, figure=None): + """draw a solution trajectory (for overlay on a contour plot) + +x is a sequence of values for one parameter (i.e. a parameter trajectory) +y is a sequence of values for one parameter (i.e. a parameter trajectory) +cost is the solution trajectory (i.e. costs); if provided, plot a 3D contour +if scale is provided, scale the intensity as 'z = log(4*z*scale+1)+2' +if shift is provided, shift the intensity as 'z = z+shift' (useful for -z's) +if style is provided, set the line style (e.g. 'w-o', 'k-', 'ro') +if figure is provided, plot to an existing figure + """ + if not figure: figure = plt.figure() + + if cost: kwds = {'projection':'3d'} # 3D + else: kwds = {} # 2D + ax = figure.gca(**kwds) + + if style in [None, False]: + style = 'w-o' #if not scale else 'k-o' + if cost: # is 3D, cost is needed + import numpy + if shift: + if shift is True: #NOTE: MAY NOT be the exact minimum + shift = max(-numpy.min(cost), 0.0) + 0.5 # a good guess + cost = numpy.asarray(cost)+shift + if scale: + cost = numpy.asarray(cost) + cost = numpy.log(4*cost*scale+1)+2 + ax.plot(x,y,cost, style, linewidth=2, markersize=4) + #XXX: need to 'correct' the z-axis (or provide easy conversion) + else: # is 2D, cost not needed + ax.plot(x,y, style, linewidth=2, markersize=4) + return figure + + +def draw_slice(f, x, y=None, scale=True, shift=False): + """plot a slice of a 2D function 'f' in 1D + +x is an array used to set up the axis +y is a fixed value for the 2nd axis +if scale is provided, scale the intensity as 'z = log(4*z*scale+1)+2' +if shift is provided, shift the intensity as 'z = z+shift' (useful for -z's) + +NOTE: when plotting the 'y-axis' at fixed 'x', +pass the array to 'y' and the fixed value to 'x' + """ + import numpy + + if y is None: + y = 0.0 + x, y = numpy.meshgrid(x, y) + plotx = True if numpy.all(y == y[0,0]) else False + + z = 0*x + s,t = x.shape + for i in range(s): + for j in range(t): + xx,yy = x[i,j], y[i,j] + z[i,j] = f([xx,yy]) + if shift: + if shift is True: shift = max(-numpy.min(z), 0.0) + 0.5 # exact minimum + z = z+shift + if scale: z = numpy.log(4*z*scale+1)+2 + #XXX: need to 'correct' the z-axis (or provide easy conversion) + + fig = plt.figure() + ax = fig.gca() + ax.autoscale(tight=True) + if plotx: + ax.plot(x.reshape(-1), z.reshape(-1)) + else: + ax.plot(y.reshape(-1), z.reshape(-1)) + return fig + + +def draw_contour(f, x, y=None, surface=False, fill=True, scale=True, shift=False, density=5): + """draw a contour plot for a given 2D function 'f' + +x and y are arrays used to set up a 2D mesh grid +if fill is True, color fill the contours +if surface is True, plot the contours as a 3D projection +if scale is provided, scale the intensity as 'z = log(4*z*scale+1)+2' +if shift is provided, shift the intensity as 'z = z+shift' (useful for -z's) +use density to adjust the number of contour lines + """ + import numpy + + if y is None: + y = x + x, y = numpy.meshgrid(x, y) + + z = 0*x + s,t = x.shape + for i in range(s): + for j in range(t): + xx,yy = x[i,j], y[i,j] + z[i,j] = f([xx,yy]) + if shift: + if shift is True: shift = max(-numpy.min(z), 0.0) + 0.5 # exact minimum + z = z+shift + if scale: z = numpy.log(4*z*scale+1)+2 + #XXX: need to 'correct' the z-axis (or provide easy conversion) + + fig = plt.figure() + if surface and fill is None: # 'hidden' option; full 3D surface plot + ax = fig.gca(projection='3d') + d = max(11 - density, 1) # or 1/density ? + kwds = {'rstride':d,'cstride':d,'cmap':cm.jet,'linewidth':0} + ax.plot_surface(x, y, z, **kwds) + else: + if surface: kwds = {'projection':'3d'} # 3D + elif surface is None: # 1D + raise NotImplementedError('need to add an option string parser') + else: kwds = {} # 2D + ax = fig.gca(**kwds) + density = 10*density + if fill: plotf = ax.contourf # filled contours + else: plotf = ax.contour # wire contours + plotf(x, y, z, density, cmap=cm.jet) + return fig + + +def model_plotter(cmdargs=[]): + """convenience function providing a function interface to the model + plotter. The ``cmdargs`` can either be a :class:`basestring` or + a :class:`list`. + + See ``mystic.model_plotter`` module documentation for a complete + description of the model plotter and its associated arguments. + + Examples: + + >>> from mystic.model_plotter import model_plotter + >>> model_plotter(cmdargs='mystic.models.zimmermann log.txt -b "-5:10:.1, -5:10:.1" -d -x 1') + + or + + >>> from mystic.model_plotter import model_plotter + >>> model_plotter(cmdargs=['mystic.models.zimmermann', 'log.txt', '-b', '-5:10:.1, -5:10:.1', '-d', '-x', '1']) + + """ + + if isinstance(cmdargs, basestring): + import shlex + cmdargs = shlex.split(cmdargs) + + if not isinstance(cmdargs, list): + raise Exception("'cmdargs' is expected to either be a list or a string") + + #FIXME: should be able to: + # - apply a constraint as a region of NaN -- apply when 'xx,yy=x[ij],y[ij]' + # - apply a penalty by shifting the surface (plot w/alpha?) -- as above + # - build an appropriately-sized default grid (from logfile info) + # - move all mulit-id param/cost reading into read_history + #FIXME: current issues: + # - 1D slice and projection work for 2D function, but aren't "pretty" + # - 1D slice and projection for 1D function, is it meaningful and correct? + # - should be able to plot from solver.genealogy (multi-monitor?) [1D,2D,3D?] + # - should be able to scale 'z-axis' instead of scaling 'z' itself + # (see https://github.com/matplotlib/matplotlib/issues/209) + # - if trajectory outside contour grid, will increase bounds + # (see support_hypercube.py for how to fix bounds) + + #XXX: note that 'argparse' is new as of python2.7 + from optparse import OptionParser + parser = OptionParser(usage=__doc__) + parser.add_option("-p","--plot-filepath",action="store",dest="plot_filepath",\ + metavar="STR",default=None, + help="save generated plot") + parser.add_option("-b","--bounds",action="store",dest="bounds",\ + metavar="STR",default="-5:5:.1, -5:5:.1", + help="indicator string to set plot bounds and density") + parser.add_option("-l","--label",action="store",dest="label",\ + metavar="STR",default=",,", + help="string to assign label to axis") + parser.add_option("-n","--nid",action="store",dest="id",\ + metavar="INT",default=None, + help="id # of the nth simultaneous points to plot") + parser.add_option("-i","--iter",action="store",dest="stop",\ + metavar="STR",default=":", + help="string for smallest:largest iterations to plot") + parser.add_option("-r","--reduce",action="store",dest="reducer",\ + metavar="STR",default="None", + help="import path of output reducer function") + parser.add_option("-x","--scale",action="store",dest="zscale",\ + metavar="INT",default=0.0, + help="scale plotted cost by z=log(4*z*scale+1)+2") + parser.add_option("-z","--shift",action="store",dest="zshift",\ + metavar="INT",default=0.0, + help="shift plotted cost by z=z+shift") + parser.add_option("-f","--fill",action="store_true",dest="fill",\ + default=False,help="plot using filled contours") + parser.add_option("-d","--depth",action="store_true",dest="surface",\ + default=False,help="plot contours showing depth in 3D") + parser.add_option("-o","--dots",action="store_true",dest="dots",\ + default=False,help="show trajectory points in plot") + parser.add_option("-j","--join",action="store_true",dest="line",\ + default=False,help="connect trajectory points in plot") + parsed_opts, parsed_args = parser.parse_args(cmdargs) + + # get the import path for the model + model = parsed_args[0] # e.g. 'mystic.models.rosen' + if "None" == model: model = None #XXX: 'required'... allow this? + + try: # get the name of the parameter log file + source = parsed_args[1] # e.g. 'log.txt' + except: + source = None + + try: # select the bounds + options = parsed_opts.bounds # format is "-1:10:.1, -1:10:.1, 1.0" + except: + options = "-5:5:.1, -5:5:.1" + + try: # plot using filled contours + fill = parsed_opts.fill + except: + fill = False + + try: # plot contours showing depth in 3D + surface = parsed_opts.surface + except: + surface = False + + #XXX: can't do '-x' with no argument given (use T/F instead?) + try: # scale plotted cost by z=log(4*z*scale+1)+2 + scale = float(parsed_opts.zscale) + if not scale: scale = False + except: + scale = False + + #XXX: can't do '-z' with no argument given + try: # shift plotted cost by z=z+shift + shift = float(parsed_opts.zshift) + if not shift: shift = False + except: + shift = False + + try: # import path of output reducer function + reducer = parsed_opts.reducer # e.g. 'numpy.add' + if "None" == reducer: reducer = None + except: + reducer = None + + style = '-' # default linestyle + if parsed_opts.dots: + mark = 'o' # marker=mark + # when using 'dots', also can turn off 'line' + if not parsed_opts.line: + style = '' # linestyle='None' + else: + mark = '' + color = 'w' if fill else 'k' + style = color + style + mark + + try: # select labels for the axes + label = parsed_opts.label.split(',') # format is "x, y, z" + except: + label = ['','',''] + + try: # select which 'id' to plot results for + ids = (int(parsed_opts.id),) #XXX: allow selecting more than one id ? + except: + ids = None # i.e. 'all' + + try: # select which iteration to stop plotting at + stop = parsed_opts.stop # format is "1:10:1" + stop = stop if ":" in stop else ":"+stop + except: + stop = ":" + + ################################################# + solver = None # set to 'mystic.solvers.fmin' (or similar) for 'live' fits + #NOTE: 'live' runs constrain params explicitly in the solver, then reduce + # dimensions appropriately so results can be 2D contour plotted. + # When working with legacy results that have more than 2 params, + # the trajectory WILL NOT follow the masked surface generated + # because the masked params were NOT fixed when the solver was run. + ################################################# + + from mystic.tools import reduced, masked, partial + + # process inputs + select, spec, mask = parse_input(options) + x,y = parse_axes(spec, grid=True) # grid=False for 1D plots + #FIXME: does grid=False still make sense here...? + if reducer: reducer = get_instance(reducer) + if solver and (not source or not model): + raise RuntimeError('a model and results filename are required') + elif not source and not model: + raise RuntimeError('a model or a results file is required') + if model: + model = get_instance(model) + # need a reducer if model returns an array + if reducer: model = reduced(reducer, arraylike=False)(model) + + if solver: + # if 'live'... pick a solver + solver = 'mystic.solvers.fmin' + solver = get_instance(solver) + xlen = len(select)+len(mask) + if solver.__name__.startswith('diffev'): + initial = [(-1,1)]*xlen + else: + initial = [0]*xlen + from mystic.monitors import VerboseLoggingMonitor + itermon = VerboseLoggingMonitor(filename=source, new=True) + # explicitly constrain parameters + model = partial(mask)(model) + # solve + sol = solver(model, x0=initial, itermon=itermon) + + #-OVERRIDE-INPUTS-# + import numpy + # read trajectories from monitor (comment out to use logfile) + source = itermon + # if negative minimum, shift by the 'solved minimum' plus an epsilon + shift = max(-numpy.min(itermon.y), 0.0) + 0.5 # a good guess + #-----------------# + + if model: # for plotting, implicitly constrain by reduction + model = masked(mask)(model) + + ## plot the surface in 1D + #if solver: v=sol[-1] + #elif source: v=cost[-1] + #else: v=None + #fig0 = draw_slice(model, x=x, y=v, scale=scale, shift=shift) + # plot the surface in 2D or 3D + fig = draw_contour(model, x, y, surface=surface, fill=fill, scale=scale, shift=shift) + else: + #fig0 = None + fig = None + + if source: + # params are the parameter trajectories + # cost is the solution trajectory + params, cost = get_history(source, ids) + if len(cost) > 1: style = style[1:] # 'auto-color' #XXX: or grayscale? + + for p,c in zip(params, cost): + ## project trajectory on a 1D slice of model surface #XXX: useful? + #s = select[0] if len(select) else 0 + #px = p[int(s)] # draw_projection requires one parameter + ## ignore everything after 'stop' + #_c = eval('c[%s]' % stop) + #_x = eval('px[%s]' % stop) + #fig0 = draw_projection(_x,_c, style=style, scale=scale, shift=shift, figure=fig0) + + # plot the trajectory on the model surface (2D or 3D) + # get two selected params #XXX: what if len(select)<2? or len(p)<2? + p = [p[int(i)] for i in select[:2]] + px,py = p # draw_trajectory requires two parameters + # ignore everything after 'stop' + _x = eval('px[%s]' % stop) + _y = eval('py[%s]' % stop) + _c = eval('c[%s]' % stop) if surface else None + fig = draw_trajectory(_x,_y,_c, style=style, scale=scale, shift=shift, figure=fig) + + # add labels to the axes + if surface: kwds = {'projection':'3d'} # 3D + else: kwds = {} # 2D + ax = fig.gca(**kwds) + ax.set_xlabel(label[0]) + ax.set_ylabel(label[1]) + if surface: ax.set_zlabel(label[2]) + + if not parsed_opts.plot_filepath: + plt.show() + else: + fig.savefig(parsed_opts.plot_filepath) + + +if __name__=='__main__': + pass + +# End of file diff --git a/scripts/mystic_model_plotter.py b/scripts/mystic_model_plotter.py index 9893376c..10d922a0 100644 --- a/scripts/mystic_model_plotter.py +++ b/scripts/mystic_model_plotter.py @@ -1,548 +1,19 @@ #!/usr/bin/env python # # Author: Mike McKerns (mmckerns @caltech and @uqfoundation) +# Author: Jean-Christophe Fillion-Robin (jchris.fillionr @kitware.com) # Copyright (c) 1997-2015 California Institute of Technology. # License: 3-clause BSD. The full license text is available at: # - http://trac.mystic.cacr.caltech.edu/project/mystic/browser/mystic/LICENSE -__doc__ = """ -mystic_model_plotter.py [options] model (filename) - -generate surface contour plots for model, specified by full import path -generate model trajectory from logfile (or solver restart file), if provided - -The option "bounds" takes an indicator string, where the bounds should -be given as comma-separated slices. For example, using bounds = "-1:10, 0:20" -will set the lower and upper bounds for x to be (-1,10) and y to be (0,20). -The "step" can also be given, to control the number of lines plotted in the -grid. Thus "-1:10:.1, 0:20" would set the bounds as above, but use increments -of .1 along x and the default step along y. For models with > 2D, the bounds -can be used to specify 2 dimensions plus fixed values for remaining dimensions. -Thus, "-1:10, 0:20, 1.0" would plot the 2D surface where the z-axis was fixed -at z=1.0. - -The option "label" takes comma-separated strings. For example, label = "x,y," -will place 'x' on the x-axis, 'y' on the y-axis, and nothing on the z-axis. -LaTeX is also accepted. For example, label = "$ h $, $ {\\alpha}$, $ v$" will -label the axes with standard LaTeX math formatting. Note that the leading -space is required, while a trailing space aligns the text with the axis -instead of the plot frame. - -The option "reduce" can be given to reduce the output of a model to a scalar, -thus converting 'model(params)' to 'reduce(model(params))'. A reducer is given -by the import path (e.g. 'numpy.add'). The option "scale" will convert the plot -to log-scale, and scale the cost by 'z=log(4*z*scale+1)+2'. This is useful for -visualizing small contour changes around the minimium. If using log-scale -produces negative numbers, the option "shift" can be used to shift the cost -by 'z=z+shift'. Both shift and scale are intended to help visualize contours. - -Required Inputs: - model full import path for the model (e.g. mystic.models.rosen) - -Additional Inputs: - filename name of the convergence logfile (e.g. log.txt) -""" - -from mpl_toolkits.mplot3d import axes3d -import matplotlib.pyplot as plt -from matplotlib import cm - -from mystic.munge import read_history -from mystic.munge import logfile_reader, raw_to_support - -#XXX: better if reads single id only? (e.g. same interface as read_history) -def get_history(source, ids=None): - """get params and cost from the given source - -source is the name of the trajectory logfile (or solver instance) -if provided, ids are the list of 'run ids' to select - """ - try: # if it's a logfile, it might be multi-id - step, param, cost = logfile_reader(source) - except: # it's not a logfile, so read and return - param, cost = read_history(source) - return [param],[cost] - - # split (i,id) into iteration and id - multinode = len(step[0]) - 1 #XXX: what if step = []? - if multinode: id = [i[1] for i in step] - else: id = [0 for i in step] - - params = [[] for i in range(max(id) + 1)] - costs = [[] for i in range(len(params))] - # populate params for each id with the corresponding (param,cost) - for i in range(len(id)): - if ids is None or id[i] in ids: # take only the selected 'id' - params[id[i]].append(param[i]) - costs[id[i]].append(cost[i]) - params = [r for r in params if len(r)] # only keep selected 'ids' - costs = [r for r in costs if len(r)] # only keep selected 'ids' - - # convert to support format - for i in range(len(params)): - params[i], costs[i] = raw_to_support(params[i], costs[i]) - return params, costs - - -def get_instance(location, *args, **kwds): - """given the import location of a model or model class, return the model - -args and kwds will be passed to the constructor of the model class - """ - package, target = location.rsplit('.',1) - exec "from %s import %s as model" % (package, target) - import inspect - if inspect.isclass(model): - model = model(*args, **kwds) - return model - - -def parse_input(option): - """parse 'option' string into 'select', 'axes', and 'mask' - -select contains the dimension specifications on which to plot -axes holds the indicies of the parameters selected to plot -mask is a dictionary of the parameter indicies and fixed values - -For example, - >>> select, axes, mask = parse_input("-1:10:.1, 0.0, 5.0, -50:50:.5") - >>> select - [0, 3] - >>> axes - "-1:10:.1, -50:50:.5" - >>> mask - {1: 0.0, 2: 5.0} - """ - option = option.split(',') - select = [] - axes = [] - mask = {} - for index,value in enumerate(option): - if ":" in value: - select.append(index) - axes.append(value) - else: - mask.update({index:float(value)}) - axes = ','.join(axes) - return select, axes, mask - - -def parse_axes(option, grid=True): - """parse option string into grid axes; using modified numpy.ogrid notation - -For example: - option='-1:10:.1, 0:10:.1' yields x,y=ogrid[-1:10:.1,0:10:.1], - -If grid is False, accept options suitable for line plotting. -For example: - option='-1:10' yields x=ogrid[-1:10] and y=0, - option='-1:10, 2' yields x=ogrid[-1:10] and y=2, - -Returns tuple (x,y) with 'x,y' defined above. - """ - import numpy - option = option.split(',') - opt = dict(zip(['x','y','z'],option)) - if len(option) > 2 or len(option) < 1: - raise ValueError("invalid format string: '%s'" % ','.join(option)) - z = bool(grid) - if len(option) == 1: opt['y'] = '0' - xd = True if ':' in opt['x'] else False - yd = True if ':' in opt['y'] else False - #XXX: accepts option='3:1', '1:1', and '1:2:10' (try to catch?) - if xd and yd: - try: # x,y form a 2D grid - exec('x,y = numpy.ogrid[%s,%s]' % (opt['x'],opt['y'])) - except: # AttributeError: - raise ValueError("invalid format string: '%s'" % ','.join(option)) - elif xd and not z: - try: - exec('x = numpy.ogrid[%s]' % opt['x']) - y = float(opt['y']) - except: # (AttributeError, SyntaxError, ValueError): - raise ValueError("invalid format string: '%s'" % ','.join(option)) - elif yd and not z: - try: - x = float(opt['x']) - exec('y = numpy.ogrid[%s]' % opt['y']) - except: # (AttributeError, SyntaxError, ValueError): - raise ValueError("invalid format string: '%s'" % ','.join(option)) - else: - raise ValueError("invalid format string: '%s'" % ','.join(option)) - if not x.size or not y.size: - raise ValueError("invalid format string: '%s'" % ','.join(option)) - return x,y - - -def draw_projection(x, cost, scale=True, shift=False, style=None, figure=None): - """draw a solution trajectory (for overlay on a 1D plot) - -x is the sequence of values for one parameter (i.e. a parameter trajectory) -cost is the sequence of costs (i.e. the solution trajectory) -if scale is provided, scale the intensity as 'z = log(4*z*scale+1)+2' -if shift is provided, shift the intensity as 'z = z+shift' (useful for -z's) -if style is provided, set the line style (e.g. 'w-o', 'k-', 'ro') -if figure is provided, plot to an existing figure - """ - if not figure: figure = plt.figure() - ax = figure.gca() - ax.autoscale(tight=True) - - if style in [None, False]: - style = 'k-o' - import numpy - if shift: - if shift is True: #NOTE: MAY NOT be the exact minimum - shift = max(-numpy.min(cost), 0.0) + 0.5 # a good guess - cost = numpy.asarray(cost)+shift - cost = numpy.asarray(cost) - if scale: - cost = numpy.log(4*cost*scale+1)+2 - - ax.plot(x,cost, style, linewidth=2, markersize=4) - #XXX: need to 'correct' the z-axis (or provide easy conversion) - return figure - - -def draw_trajectory(x, y, cost=None, scale=True, shift=False, style=None, figure=None): - """draw a solution trajectory (for overlay on a contour plot) - -x is a sequence of values for one parameter (i.e. a parameter trajectory) -y is a sequence of values for one parameter (i.e. a parameter trajectory) -cost is the solution trajectory (i.e. costs); if provided, plot a 3D contour -if scale is provided, scale the intensity as 'z = log(4*z*scale+1)+2' -if shift is provided, shift the intensity as 'z = z+shift' (useful for -z's) -if style is provided, set the line style (e.g. 'w-o', 'k-', 'ro') -if figure is provided, plot to an existing figure - """ - if not figure: figure = plt.figure() - - if cost: kwds = {'projection':'3d'} # 3D - else: kwds = {} # 2D - ax = figure.gca(**kwds) - - if style in [None, False]: - style = 'w-o' #if not scale else 'k-o' - if cost: # is 3D, cost is needed - import numpy - if shift: - if shift is True: #NOTE: MAY NOT be the exact minimum - shift = max(-numpy.min(cost), 0.0) + 0.5 # a good guess - cost = numpy.asarray(cost)+shift - if scale: - cost = numpy.asarray(cost) - cost = numpy.log(4*cost*scale+1)+2 - ax.plot(x,y,cost, style, linewidth=2, markersize=4) - #XXX: need to 'correct' the z-axis (or provide easy conversion) - else: # is 2D, cost not needed - ax.plot(x,y, style, linewidth=2, markersize=4) - return figure - - -def draw_slice(f, x, y=None, scale=True, shift=False): - """plot a slice of a 2D function 'f' in 1D - -x is an array used to set up the axis -y is a fixed value for the 2nd axis -if scale is provided, scale the intensity as 'z = log(4*z*scale+1)+2' -if shift is provided, shift the intensity as 'z = z+shift' (useful for -z's) - -NOTE: when plotting the 'y-axis' at fixed 'x', -pass the array to 'y' and the fixed value to 'x' - """ - import numpy - - if y is None: - y = 0.0 - x, y = numpy.meshgrid(x, y) - plotx = True if numpy.all(y == y[0,0]) else False - - z = 0*x - s,t = x.shape - for i in range(s): - for j in range(t): - xx,yy = x[i,j], y[i,j] - z[i,j] = f([xx,yy]) - if shift: - if shift is True: shift = max(-numpy.min(z), 0.0) + 0.5 # exact minimum - z = z+shift - if scale: z = numpy.log(4*z*scale+1)+2 - #XXX: need to 'correct' the z-axis (or provide easy conversion) - - fig = plt.figure() - ax = fig.gca() - ax.autoscale(tight=True) - if plotx: - ax.plot(x.reshape(-1), z.reshape(-1)) - else: - ax.plot(y.reshape(-1), z.reshape(-1)) - return fig - - -def draw_contour(f, x, y=None, surface=False, fill=True, scale=True, shift=False, density=5): - """draw a contour plot for a given 2D function 'f' - -x and y are arrays used to set up a 2D mesh grid -if fill is True, color fill the contours -if surface is True, plot the contours as a 3D projection -if scale is provided, scale the intensity as 'z = log(4*z*scale+1)+2' -if shift is provided, shift the intensity as 'z = z+shift' (useful for -z's) -use density to adjust the number of contour lines - """ - import numpy - - if y is None: - y = x - x, y = numpy.meshgrid(x, y) - - z = 0*x - s,t = x.shape - for i in range(s): - for j in range(t): - xx,yy = x[i,j], y[i,j] - z[i,j] = f([xx,yy]) - if shift: - if shift is True: shift = max(-numpy.min(z), 0.0) + 0.5 # exact minimum - z = z+shift - if scale: z = numpy.log(4*z*scale+1)+2 - #XXX: need to 'correct' the z-axis (or provide easy conversion) - - fig = plt.figure() - if surface and fill is None: # 'hidden' option; full 3D surface plot - ax = fig.gca(projection='3d') - d = max(11 - density, 1) # or 1/density ? - kwds = {'rstride':d,'cstride':d,'cmap':cm.jet,'linewidth':0} - ax.plot_surface(x, y, z, **kwds) - else: - if surface: kwds = {'projection':'3d'} # 3D - elif surface is None: # 1D - raise NotImplementedError('need to add an option string parser') - else: kwds = {} # 2D - ax = fig.gca(**kwds) - density = 10*density - if fill: plotf = ax.contourf # filled contours - else: plotf = ax.contour # wire contours - plotf(x, y, z, density, cmap=cm.jet) - return fig +from mystic.model_plotter import model_plotter +__doc__ = model_plotter.__doc__ if __name__ == '__main__': - #FIXME: should be able to: - # - apply a constraint as a region of NaN -- apply when 'xx,yy=x[ij],y[ij]' - # - apply a penalty by shifting the surface (plot w/alpha?) -- as above - # - build an appropriately-sized default grid (from logfile info) - # - move all mulit-id param/cost reading into read_history - #FIXME: current issues: - # - 1D slice and projection work for 2D function, but aren't "pretty" - # - 1D slice and projection for 1D function, is it meaningful and correct? - # - should be able to plot from solver.genealogy (multi-monitor?) [1D,2D,3D?] - # - should be able to scale 'z-axis' instead of scaling 'z' itself - # (see https://github.com/matplotlib/matplotlib/issues/209) - # - if trajectory outside contour grid, will increase bounds - # (see support_hypercube.py for how to fix bounds) - - #XXX: note that 'argparse' is new as of python2.7 - from optparse import OptionParser - parser = OptionParser(usage=__doc__) - parser.add_option("-b","--bounds",action="store",dest="bounds",\ - metavar="STR",default="-5:5:.1, -5:5:.1", - help="indicator string to set plot bounds and density") - parser.add_option("-l","--label",action="store",dest="label",\ - metavar="STR",default=",,", - help="string to assign label to axis") - parser.add_option("-n","--nid",action="store",dest="id",\ - metavar="INT",default=None, - help="id # of the nth simultaneous points to plot") - parser.add_option("-i","--iter",action="store",dest="stop",\ - metavar="STR",default=":", - help="string for smallest:largest iterations to plot") - parser.add_option("-r","--reduce",action="store",dest="reducer",\ - metavar="STR",default="None", - help="import path of output reducer function") - parser.add_option("-x","--scale",action="store",dest="zscale",\ - metavar="INT",default=0.0, - help="scale plotted cost by z=log(4*z*scale+1)+2") - parser.add_option("-z","--shift",action="store",dest="zshift",\ - metavar="INT",default=0.0, - help="shift plotted cost by z=z+shift") - parser.add_option("-f","--fill",action="store_true",dest="fill",\ - default=False,help="plot using filled contours") - parser.add_option("-d","--depth",action="store_true",dest="surface",\ - default=False,help="plot contours showing depth in 3D") - parser.add_option("-o","--dots",action="store_true",dest="dots",\ - default=False,help="show trajectory points in plot") - parser.add_option("-j","--join",action="store_true",dest="line",\ - default=False,help="connect trajectory points in plot") - parsed_opts, parsed_args = parser.parse_args() - - # get the import path for the model - model = parsed_args[0] # e.g. 'mystic.models.rosen' - if "None" == model: model = None #XXX: 'required'... allow this? - - try: # get the name of the parameter log file - source = parsed_args[1] # e.g. 'log.txt' - except: - source = None - - try: # select the bounds - options = parsed_opts.bounds # format is "-1:10:.1, -1:10:.1, 1.0" - except: - options = "-5:5:.1, -5:5:.1" - - try: # plot using filled contours - fill = parsed_opts.fill - except: - fill = False - - try: # plot contours showing depth in 3D - surface = parsed_opts.surface - except: - surface = False - - #XXX: can't do '-x' with no argument given (use T/F instead?) - try: # scale plotted cost by z=log(4*z*scale+1)+2 - scale = float(parsed_opts.zscale) - if not scale: scale = False - except: - scale = False - - #XXX: can't do '-z' with no argument given - try: # shift plotted cost by z=z+shift - shift = float(parsed_opts.zshift) - if not shift: shift = False - except: - shift = False - - try: # import path of output reducer function - reducer = parsed_opts.reducer # e.g. 'numpy.add' - if "None" == reducer: reducer = None - except: - reducer = None - - style = '-' # default linestyle - if parsed_opts.dots: - mark = 'o' # marker=mark - # when using 'dots', also can turn off 'line' - if not parsed_opts.line: - style = '' # linestyle='None' - else: - mark = '' - color = 'w' if fill else 'k' - style = color + style + mark - - try: # select labels for the axes - label = parsed_opts.label.split(',') # format is "x, y, z" - except: - label = ['','',''] - - try: # select which 'id' to plot results for - ids = (int(parsed_opts.id),) #XXX: allow selecting more than one id ? - except: - ids = None # i.e. 'all' - - try: # select which iteration to stop plotting at - stop = parsed_opts.stop # format is "1:10:1" - stop = stop if ":" in stop else ":"+stop - except: - stop = ":" - - ################################################# - solver = None # set to 'mystic.solvers.fmin' (or similar) for 'live' fits - #NOTE: 'live' runs constrain params explicitly in the solver, then reduce - # dimensions appropriately so results can be 2D contour plotted. - # When working with legacy results that have more than 2 params, - # the trajectory WILL NOT follow the masked surface generated - # because the masked params were NOT fixed when the solver was run. - ################################################# - - from mystic.tools import reduced, masked, partial - - # process inputs - select, spec, mask = parse_input(options) - x,y = parse_axes(spec, grid=True) # grid=False for 1D plots - #FIXME: does grid=False still make sense here...? - if reducer: reducer = get_instance(reducer) - if solver and (not source or not model): - raise RuntimeError('a model and results filename are required') - elif not source and not model: - raise RuntimeError('a model or a results file is required') - if model: - model = get_instance(model) - # need a reducer if model returns an array - if reducer: model = reduced(reducer, arraylike=False)(model) - - if solver: - # if 'live'... pick a solver - solver = 'mystic.solvers.fmin' - solver = get_instance(solver) - xlen = len(select)+len(mask) - if solver.__name__.startswith('diffev'): - initial = [(-1,1)]*xlen - else: - initial = [0]*xlen - from mystic.monitors import VerboseLoggingMonitor - itermon = VerboseLoggingMonitor(filename=source, new=True) - # explicitly constrain parameters - model = partial(mask)(model) - # solve - sol = solver(model, x0=initial, itermon=itermon) - - #-OVERRIDE-INPUTS-# - import numpy - # read trajectories from monitor (comment out to use logfile) - source = itermon - # if negative minimum, shift by the 'solved minimum' plus an epsilon - shift = max(-numpy.min(itermon.y), 0.0) + 0.5 # a good guess - #-----------------# - - if model: # for plotting, implicitly constrain by reduction - model = masked(mask)(model) - - ## plot the surface in 1D - #if solver: v=sol[-1] - #elif source: v=cost[-1] - #else: v=None - #fig0 = draw_slice(model, x=x, y=v, scale=scale, shift=shift) - # plot the surface in 2D or 3D - fig = draw_contour(model, x, y, surface=surface, fill=fill, scale=scale, shift=shift) - else: - #fig0 = None - fig = None - - if source: - # params are the parameter trajectories - # cost is the solution trajectory - params, cost = get_history(source, ids) - if len(cost) > 1: style = style[1:] # 'auto-color' #XXX: or grayscale? - - for p,c in zip(params, cost): - ## project trajectory on a 1D slice of model surface #XXX: useful? - #s = select[0] if len(select) else 0 - #px = p[int(s)] # draw_projection requires one parameter - ## ignore everything after 'stop' - #_c = eval('c[%s]' % stop) - #_x = eval('px[%s]' % stop) - #fig0 = draw_projection(_x,_c, style=style, scale=scale, shift=shift, figure=fig0) - - # plot the trajectory on the model surface (2D or 3D) - # get two selected params #XXX: what if len(select)<2? or len(p)<2? - p = [p[int(i)] for i in select[:2]] - px,py = p # draw_trajectory requires two parameters - # ignore everything after 'stop' - _x = eval('px[%s]' % stop) - _y = eval('py[%s]' % stop) - _c = eval('c[%s]' % stop) if surface else None - fig = draw_trajectory(_x,_y,_c, style=style, scale=scale, shift=shift, figure=fig) - - # add labels to the axes - if surface: kwds = {'projection':'3d'} # 3D - else: kwds = {} # 2D - ax = fig.gca(**kwds) - ax.set_xlabel(label[0]) - ax.set_ylabel(label[1]) - if surface: ax.set_zlabel(label[2]) - plt.show() + import sys + model_plotter(cmdargs=sys.argv[1:]) # EOF