New files

GPy/plotting/matplot_dep/__init__.py

@@ -0,0 +1,16 @@
+# Copyright (c) 2014, GPy authors (see AUTHORS.txt).
+# Licensed under the BSD 3-clause license (see LICENSE.txt)
+
+import base_plots
+import models_plots
+import priors_plots
+import variational_plots
+import kernel_plots
+import svig_plots
+import dim_reduction_plots
+import mapping_plots
+import Tango
+import visualize
+import latent_space_visualizations
+import netpbmfile
+import inference_plots

GPy/plotting/matplot_dep/inference_plots.py

@@ -0,0 +1,28 @@
+# Copyright (c) 2012, GPy authors (see AUTHORS.txt).
+# Licensed under the BSD 3-clause license (see LICENSE.txt)
+
+import pylab as pb
+#import numpy as np
+#import Tango
+#from base_plots import gpplot, x_frame1D, x_frame2D
+
+
+def plot_optimizer(optimizer):
+    if optimizer.trace == None:
+        print "No trace present so I can't plot it. Please check that the optimizer actually supplies a trace."
+    else:
+        pb.figure()
+        pb.plot(optimizer.trace)
+        pb.xlabel('Iteration')
+        pb.ylabel('f(x)')
+
+def plot_sgd_traces(optimizer):
+    pb.figure()
+    pb.subplot(211)
+    pb.title('Parameters')
+    for k in optimizer.param_traces.keys():
+        pb.plot(optimizer.param_traces[k], label=k)
+    pb.legend(loc=0)
+    pb.subplot(212)
+    pb.title('Objective function')
+    pb.plot(optimizer.fopt_trace)

GPy/plotting/matplot_dep/kernel_plots.py

@@ -0,0 +1,137 @@
+# Copyright (c) 2012, GPy authors (see AUTHORS.txt).
+# Licensed under the BSD 3-clause license (see LICENSE.txt)
+
+import sys
+import numpy as np
+import pylab as pb
+import Tango
+from matplotlib.textpath import TextPath
+from matplotlib.transforms import offset_copy
+
+
+def plot_ARD(kernel, fignum=None, ax=None, title='', legend=False):
+    """If an ARD kernel is present, plot a bar representation using matplotlib
+
+    :param fignum: figure number of the plot
+    :param ax: matplotlib axis to plot on
+    :param title:
+        title of the plot,
+        pass '' to not print a title
+        pass None for a generic title
+    """
+    if ax is None:
+        fig = pb.figure(fignum)
+        ax = fig.add_subplot(111)
+    else:
+        fig = ax.figure
+    Tango.reset()
+    xticklabels = []
+    bars = []
+    x0 = 0
+    for p in kernel._parameters_:
+        c = Tango.nextMedium()
+        if hasattr(p, 'ARD') and p.ARD:
+            if title is None:
+                ax.set_title('ARD parameters, %s kernel' % p.name)
+            else:
+                ax.set_title(title)
+            if isinstance(p, Linear):
+                ard_params = p.variances
+            else:
+                ard_params = 1. / p.lengthscale
+
+            x = np.arange(x0, x0 + len(ard_params))
+            bars.append(ax.bar(x, ard_params, align='center', color=c, edgecolor='k', linewidth=1.2, label=p.name.replace("_"," ")))
+            xticklabels.extend([r"$\mathrm{{{name}}}\ {x}$".format(name=p.name, x=i) for i in np.arange(len(ard_params))])
+            x0 += len(ard_params)
+    x = np.arange(x0)
+    transOffset = offset_copy(ax.transData, fig=fig,
+                              x=0., y= -2., units='points')
+    transOffsetUp = offset_copy(ax.transData, fig=fig,
+                              x=0., y=1., units='points')
+    for bar in bars:
+        for patch, num in zip(bar.patches, np.arange(len(bar.patches))):
+            height = patch.get_height()
+            xi = patch.get_x() + patch.get_width() / 2.
+            va = 'top'
+            c = 'w'
+            t = TextPath((0, 0), "${xi}$".format(xi=xi), rotation=0, usetex=True, ha='center')
+            transform = transOffset
+            if patch.get_extents().height <= t.get_extents().height + 3:
+                va = 'bottom'
+                c = 'k'
+                transform = transOffsetUp
+            ax.text(xi, height, "${xi}$".format(xi=int(num)), color=c, rotation=0, ha='center', va=va, transform=transform)
+    # for xi, t in zip(x, xticklabels):
+    #    ax.text(xi, maxi / 2, t, rotation=90, ha='center', va='center')
+    # ax.set_xticklabels(xticklabels, rotation=17)
+    ax.set_xticks([])
+    ax.set_xlim(-.5, x0 - .5)
+    if legend:
+        if title is '':
+            mode = 'expand'
+            if len(bars) > 1:
+                mode = 'expand'
+            ax.legend(bbox_to_anchor=(0., 1.02, 1., 1.02), loc=3,
+                      ncol=len(bars), mode=mode, borderaxespad=0.)
+            fig.tight_layout(rect=(0, 0, 1, .9))
+        else:
+            ax.legend()
+    return ax
+
+
+def plot(kernel, x=None, plot_limits=None, which_parts='all', resolution=None, *args, **kwargs):
+    if which_parts == 'all':
+        which_parts = [True] * kernel.size
+    if kernel.input_dim == 1:
+        if x is None:
+            x = np.zeros((1, 1))
+        else:
+            x = np.asarray(x)
+            assert x.size == 1, "The size of the fixed variable x is not 1"
+            x = x.reshape((1, 1))
+
+        if plot_limits == None:
+            xmin, xmax = (x - 5).flatten(), (x + 5).flatten()
+        elif len(plot_limits) == 2:
+            xmin, xmax = plot_limits
+        else:
+            raise ValueError, "Bad limits for plotting"
+
+        Xnew = np.linspace(xmin, xmax, resolution or 201)[:, None]
+        Kx = kernel.K(Xnew, x, which_parts)
+        pb.plot(Xnew, Kx, *args, **kwargs)
+        pb.xlim(xmin, xmax)
+        pb.xlabel("x")
+        pb.ylabel("k(x,%0.1f)" % x)
+
+    elif kernel.input_dim == 2:
+        if x is None:
+            x = np.zeros((1, 2))
+        else:
+            x = np.asarray(x)
+            assert x.size == 2, "The size of the fixed variable x is not 2"
+            x = x.reshape((1, 2))
+
+        if plot_limits == None:
+            xmin, xmax = (x - 5).flatten(), (x + 5).flatten()
+        elif len(plot_limits) == 2:
+            xmin, xmax = plot_limits
+        else:
+            raise ValueError, "Bad limits for plotting"
+
+        resolution = resolution or 51
+        xx, yy = np.mgrid[xmin[0]:xmax[0]:1j * resolution, xmin[1]:xmax[1]:1j * resolution]
+        xg = np.linspace(xmin[0], xmax[0], resolution)
+        yg = np.linspace(xmin[1], xmax[1], resolution)
+        Xnew = np.vstack((xx.flatten(), yy.flatten())).T
+        Kx = kernel.K(Xnew, x, which_parts)
+        Kx = Kx.reshape(resolution, resolution).T
+        pb.contour(xg, yg, Kx, vmin=Kx.min(), vmax=Kx.max(), cmap=pb.cm.jet, *args, **kwargs) # @UndefinedVariable
+        pb.xlim(xmin[0], xmax[0])
+        pb.ylim(xmin[1], xmax[1])
+        pb.xlabel("x1")
+        pb.ylabel("x2")
+        pb.title("k(x1,x2 ; %0.1f,%0.1f)" % (x[0, 0], x[0, 1]))
+    else:
+        raise NotImplementedError, "Cannot plot a kernel with more than two input dimensions"

GPy/plotting/matplot_dep/mapping_plots.py

@@ -0,0 +1,81 @@
+# Copyright (c) 2012, GPy authors (see AUTHORS.txt).
+# Licensed under the BSD 3-clause license (see LICENSE.txt)
+
+import pylab as pb
+import numpy as np
+import Tango
+from base_plots import x_frame1D, x_frame2D
+
+
+def plot_mapping(self, plot_limits=None, which_data='all', which_parts='all', resolution=None, levels=20, samples=0, fignum=None, ax=None, fixed_inputs=[], linecol=Tango.colorsHex['darkBlue']):
+    """
+    Plots the mapping associated with the model.
+      - In one dimension, the function is plotted.
+      - In two dimsensions, a contour-plot shows the function
+      - In higher dimensions, we've not implemented this yet !TODO!
+
+    Can plot only part of the data and part of the posterior functions
+    using which_data and which_functions
+
+    :param plot_limits: The limits of the plot. If 1D [xmin,xmax], if 2D [[xmin,ymin],[xmax,ymax]]. Defaluts to data limits
+    :type plot_limits: np.array
+    :param which_data: which if the training data to plot (default all)
+    :type which_data: 'all' or a slice object to slice self.X, self.Y
+    :param which_parts: which of the kernel functions to plot (additively)
+    :type which_parts: 'all', or list of bools
+    :param resolution: the number of intervals to sample the GP on. Defaults to 200 in 1D and 50 (a 50x50 grid) in 2D
+    :type resolution: int
+    :param levels: number of levels to plot in a contour plot.
+    :type levels: int
+    :param samples: the number of a posteriori samples to plot
+    :type samples: int
+    :param fignum: figure to plot on.
+    :type fignum: figure number
+    :param ax: axes to plot on.
+    :type ax: axes handle
+    :param fixed_inputs: a list of tuple [(i,v), (i,v)...], specifying that input index i should be set to value v.
+    :type fixed_inputs: a list of tuples
+    :param linecol: color of line to plot.
+    :type linecol:
+    :param levels: for 2D plotting, the number of contour levels to use is ax is None, create a new figure
+
+    """
+    # TODO include samples
+    if which_data == 'all':
+        which_data = slice(None)
+
+    if ax is None:
+        fig = pb.figure(num=fignum)
+        ax = fig.add_subplot(111)
+
+    plotdims = self.input_dim - len(fixed_inputs)
+
+    if plotdims == 1:
+
+        Xu = self.X * self._Xscale + self._Xoffset # NOTE self.X are the normalized values now
+
+        fixed_dims = np.array([i for i,v in fixed_inputs])
+        freedim = np.setdiff1d(np.arange(self.input_dim),fixed_dims)
+
+        Xnew, xmin, xmax = x_frame1D(Xu[:,freedim], plot_limits=plot_limits)
+        Xgrid = np.empty((Xnew.shape[0],self.input_dim))
+        Xgrid[:,freedim] = Xnew
+        for i,v in fixed_inputs:
+            Xgrid[:,i] = v
+
+        f = self.predict(Xgrid, which_parts=which_parts)
+        for d in range(y.shape[1]):
+            ax.plot(Xnew, f[:, d], edgecol=linecol)
+
+    elif self.X.shape[1] == 2:
+        resolution = resolution or 50
+        Xnew, _, _, xmin, xmax = x_frame2D(self.X, plot_limits, resolution)
+        x, y = np.linspace(xmin[0], xmax[0], resolution), np.linspace(xmin[1], xmax[1], resolution)
+        f = self.predict(Xnew, which_parts=which_parts)
+        m = m.reshape(resolution, resolution).T
+        ax.contour(x, y, f, levels, vmin=m.min(), vmax=m.max(), cmap=pb.cm.jet) # @UndefinedVariable
+        ax.set_xlim(xmin[0], xmax[0])
+        ax.set_ylim(xmin[1], xmax[1])
+
+    else:
+        raise NotImplementedError, "Cannot define a frame with more than two input dimensions"

GPy/plotting/matplot_dep/models_plots.py

@@ -0,0 +1,161 @@
+# Copyright (c) 2012, GPy authors (see AUTHORS.txt).
+# Licensed under the BSD 3-clause license (see LICENSE.txt)
+
+import pylab as pb
+import numpy as np
+import Tango
+from base_plots import gpplot, x_frame1D, x_frame2D
+
+
+def plot_fit(model, plot_limits=None, which_data_rows='all',
+        which_data_ycols='all', which_parts='all', fixed_inputs=[],
+        levels=20, samples=0, fignum=None, ax=None, resolution=None,
+        plot_raw=False,
+        linecol=Tango.colorsHex['darkBlue'],fillcol=Tango.colorsHex['lightBlue']):
+    """
+    Plot the posterior of the GP.
+      - In one dimension, the function is plotted with a shaded region identifying two standard deviations.
+      - In two dimsensions, a contour-plot shows the mean predicted function
+      - In higher dimensions, use fixed_inputs to plot the GP  with some of the inputs fixed.
+
+    Can plot only part of the data and part of the posterior functions
+    using which_data_rowsm which_data_ycols and which_parts
+
+    :param plot_limits: The limits of the plot. If 1D [xmin,xmax], if 2D [[xmin,ymin],[xmax,ymax]]. Defaluts to data limits
+    :type plot_limits: np.array
+    :param which_data_rows: which of the training data to plot (default all)
+    :type which_data_rows: 'all' or a slice object to slice model.X, model.Y
+    :param which_data_ycols: when the data has several columns (independant outputs), only plot these
+    :type which_data_rows: 'all' or a list of integers
+    :param which_parts: which of the kernel functions to plot (additively)
+    :type which_parts: 'all', or list of bools
+    :param fixed_inputs: a list of tuple [(i,v), (i,v)...], specifying that input index i should be set to value v.
+    :type fixed_inputs: a list of tuples
+    :param resolution: the number of intervals to sample the GP on. Defaults to 200 in 1D and 50 (a 50x50 grid) in 2D
+    :type resolution: int
+    :param levels: number of levels to plot in a contour plot.
+    :type levels: int
+    :param samples: the number of a posteriori samples to plot
+    :type samples: int
+    :param fignum: figure to plot on.
+    :type fignum: figure number
+    :param ax: axes to plot on.
+    :type ax: axes handle
+    :type output: integer (first output is 0)
+    :param linecol: color of line to plot.
+    :type linecol:
+    :param fillcol: color of fill
+    :param levels: for 2D plotting, the number of contour levels to use is ax is None, create a new figure
+    """
+    #deal with optional arguments
+    if which_data_rows == 'all':
+        which_data_rows = slice(None)
+    if which_data_ycols == 'all':
+        which_data_ycols = np.arange(model.output_dim)
+    if len(which_data_ycols)==0:
+        raise ValueError('No data selected for plotting')
+    if ax is None:
+        fig = pb.figure(num=fignum)
+        ax = fig.add_subplot(111)
+
+    #work out what the inputs are for plotting (1D or 2D)
+    fixed_dims = np.array([i for i,v in fixed_inputs])
+    free_dims = np.setdiff1d(np.arange(model.input_dim),fixed_dims)
+
+    #one dimensional plotting
+    if len(free_dims) == 1:
+
+        #define the frame on which to plot
+        resolution = resolution or 200
+        Xnew, xmin, xmax = x_frame1D(model.X[:,free_dims], plot_limits=plot_limits)
+        Xgrid = np.empty((Xnew.shape[0],model.input_dim))
+        Xgrid[:,free_dims] = Xnew
+        for i,v in fixed_inputs:
+            Xgrid[:,i] = v
+
+        #make a prediction on the frame and plot it
+        if plot_raw:
+            m, v = model._raw_predict(Xgrid, which_parts=which_parts)
+            lower = m - 2*np.sqrt(v)
+            upper = m + 2*np.sqrt(v)
+            Y = model.Y
+        else:
+            m, v, lower, upper = model.predict(Xgrid, which_parts=which_parts)
+            Y = model.Y
+        for d in which_data_ycols:
+            gpplot(Xnew, m[:, d], lower[:, d], upper[:, d], axes=ax, edgecol=linecol, fillcol=fillcol)
+            ax.plot(model.X[which_data_rows,free_dims], Y[which_data_rows, d], 'kx', mew=1.5)
+
+        #optionally plot some samples
+        if samples: #NOTE not tested with fixed_inputs
+            Ysim = model.posterior_samples(Xgrid, samples, which_parts=which_parts)
+            for yi in Ysim.T:
+                ax.plot(Xnew, yi[:,None], Tango.colorsHex['darkBlue'], linewidth=0.25)
+                #ax.plot(Xnew, yi[:,None], marker='x', linestyle='--',color=Tango.colorsHex['darkBlue']) #TODO apply this line for discrete outputs.
+
+        #add inducing inputs (if a sparse model is used)
+        if hasattr(model,"Z"):
+            Zu = model.Z[:,free_dims] * model._Xscale[:,free_dims] + model._Xoffset[:,free_dims]
+            ax.plot(Zu, np.zeros_like(Zu) + ax.get_ylim()[0], 'r|', mew=1.5, markersize=12)
+
+        #add error bars for uncertain (if input uncertainty is being modelled)
+        if hasattr(model,"has_uncertain_inputs"):
+            ax.errorbar(model.X[which_data, free_dims], model.likelihood.data[which_data, 0],
+                        xerr=2 * np.sqrt(model.X_variance[which_data, free_dims]),
+                        ecolor='k', fmt=None, elinewidth=.5, alpha=.5)
+
+
+        #set the limits of the plot to some sensible values
+        ymin, ymax = min(np.append(Y[which_data_rows, which_data_ycols].flatten(), lower)), max(np.append(Y[which_data_rows, which_data_ycols].flatten(), upper))
+        ymin, ymax = ymin - 0.1 * (ymax - ymin), ymax + 0.1 * (ymax - ymin)
+        ax.set_xlim(xmin, xmax)
+        ax.set_ylim(ymin, ymax)
+
+    #2D plotting
+    elif len(free_dims) == 2:
+
+        #define the frame for plotting on
+        resolution = resolution or 50
+        Xnew, _, _, xmin, xmax = x_frame2D(model.X[:,free_dims], plot_limits, resolution)
+        Xgrid = np.empty((Xnew.shape[0],model.input_dim))
+        Xgrid[:,free_dims] = Xnew
+        for i,v in fixed_inputs:
+            Xgrid[:,i] = v
+        x, y = np.linspace(xmin[0], xmax[0], resolution), np.linspace(xmin[1], xmax[1], resolution)
+
+        #predict on the frame and plot
+        if plot_raw:
+            m, _ = model._raw_predict(Xgrid, which_parts=which_parts)
+            Y = model.likelihood.Y
+        else:
+            m, _, _, _ = model.predict(Xgrid, which_parts=which_parts,sampling=False)
+            Y = model.likelihood.data
+        for d in which_data_ycols:
+            m_d = m[:,d].reshape(resolution, resolution).T
+            ax.contour(x, y, m_d, levels, vmin=m.min(), vmax=m.max(), cmap=pb.cm.jet)
+            ax.scatter(model.X[which_data_rows, free_dims[0]], model.X[which_data_rows, free_dims[1]], 40, Y[which_data_rows, d], cmap=pb.cm.jet, vmin=m.min(), vmax=m.max(), linewidth=0.)
+
+        #set the limits of the plot to some sensible values
+        ax.set_xlim(xmin[0], xmax[0])
+        ax.set_ylim(xmin[1], xmax[1])
+
+        if samples:
+            warnings.warn("Samples are rather difficult to plot for 2D inputs...")
+
+        #add inducing inputs (if a sparse model is used)
+        if hasattr(model,"Z"):
+            Zu = model.Z[:,free_dims] * model._Xscale[:,free_dims] + model._Xoffset[:,free_dims]
+            ax.plot(Zu[:,free_dims[0]], Zu[:,free_dims[1]], 'wo')
+
+    else:
+        raise NotImplementedError, "Cannot define a frame with more than two input dimensions"
+
+
+def plot_f_fit(model, *args, **kwargs):
+    """
+    Plot the GP's view of the world, where the data is normalized and before applying a likelihood.
+
+    All args and kwargs are passed on to models_plots.plot.
+    """
+    kwargs['plot_raw'] = True
+    plot(model,*args, **kwargs)

GPy/plotting/matplot_dep/priors_plots.py

@@ -0,0 +1,29 @@
+# Copyright (c) 2012, GPy authors (see AUTHORS.txt).
+# Licensed under the BSD 3-clause license (see LICENSE.txt)
+
+
+import numpy as np
+import pylab as pb
+
+
+def univariate_plot(prior):
+    rvs = prior.rvs(1000)
+    pb.hist(rvs, 100, normed=True)
+    xmin, xmax = pb.xlim()
+    xx = np.linspace(xmin, xmax, 1000)
+    pb.plot(xx, prior.pdf(xx), 'r', linewidth=2)
+
+def plot(prior):
+
+    if prior.input_dim == 2:
+        rvs = prior.rvs(200)
+        pb.plot(rvs[:, 0], rvs[:, 1], 'kx', mew=1.5)
+        xmin, xmax = pb.xlim()
+        ymin, ymax = pb.ylim()
+        xx, yy = np.mgrid[xmin:xmax:100j, ymin:ymax:100j]
+        xflat = np.vstack((xx.flatten(), yy.flatten())).T
+        zz = prior.pdf(xflat).reshape(100, 100)
+        pb.contour(xx, yy, zz, linewidths=2)
+
+    else:
+        raise NotImplementedError, "Cannot define a frame with more than two input dimensions"

GPy/plotting/matplot_dep/svig_plots.py

@@ -0,0 +1,43 @@
+# Copyright (c) 2012, James Hensman and Nicolo' Fusi
+# Licensed under the BSD 3-clause license (see LICENSE.txt)
+
+import numpy as np
+import pylab as pb
+
+
+def plot(model, ax=None, fignum=None, Z_height=None, **kwargs):
+
+    if ax is None:
+        fig = pb.figure(num=fignum)
+        ax = fig.add_subplot(111)
+
+    #horrible hack here:
+    data = model.likelihood.data.copy()
+    model.likelihood.data = model.Y
+    GP.plot(model, ax=ax, **kwargs)
+    model.likelihood.data = data
+
+    Zu = model.Z * model._Xscale + model._Xoffset
+    if model.input_dim==1:
+        ax.plot(model.X_batch, model.likelihood.data, 'gx',mew=2)
+        if Z_height is None:
+            Z_height = ax.get_ylim()[0]
+        ax.plot(Zu, np.zeros_like(Zu) + Z_height, 'r|', mew=1.5, markersize=12)
+
+    if model.input_dim==2:
+        ax.scatter(model.X[:,0], model.X[:,1], 20., model.Y[:,0], linewidth=0, cmap=pb.cm.jet)  # @UndefinedVariable
+        ax.plot(Zu[:,0], Zu[:,1], 'w^')
+
+def plot_traces(model):
+
+    pb.figure()
+    t = np.array(model._param_trace)
+    pb.subplot(2,1,1)
+    for l,ti in zip(model._get_param_names(),t.T):
+        if not l[:3]=='iip':
+            pb.plot(ti,label=l)
+    pb.legend(loc=0)
+
+    pb.subplot(2,1,2)
+    pb.plot(np.asarray(model._ll_trace),label='stochastic likelihood')
+    pb.legend(loc=0)

GPy/plotting/matplot_dep/variational_plots.py

@@ -0,0 +1,45 @@
+import pylab as pb
+
+def plot(parameterized, fignum=None, ax=None, colors=None):
+    """
+    Plot latent space X in 1D:
+
+        - if fig is given, create input_dim subplots in fig and plot in these
+        - if ax is given plot input_dim 1D latent space plots of X into each `axis`
+        - if neither fig nor ax is given create a figure with fignum and plot in there
+
+    colors:
+        colors of different latent space dimensions input_dim
+
+    """
+    if ax is None:
+        fig = pb.figure(num=fignum, figsize=(8, min(12, (2 * parameterized.means.shape[1]))))
+    if colors is None:
+        colors = pb.gca()._get_lines.color_cycle
+        pb.clf()
+    else:
+        colors = iter(colors)
+    plots = []
+    means, variances = param_to_array(parameterized.means, parameterized.variances)
+    x = np.arange(means.shape[0])
+    for i in range(means.shape[1]):
+        if ax is None:
+            a = fig.add_subplot(means.shape[1], 1, i + 1)
+        elif isinstance(ax, (tuple, list)):
+            a = ax[i]
+        else:
+            raise ValueError("Need one ax per latent dimension input_dim")
+        a.plot(means, c='k', alpha=.3)
+        plots.extend(a.plot(x, means.T[i], c=colors.next(), label=r"$\mathbf{{X_{{{}}}}}$".format(i)))
+        a.fill_between(x,
+                        means.T[i] - 2 * np.sqrt(variances.T[i]),
+                        means.T[i] + 2 * np.sqrt(variances.T[i]),
+                        facecolor=plots[-1].get_color(),
+                        alpha=.3)
+        a.legend(borderaxespad=0.)
+        a.set_xlim(x.min(), x.max())
+        if i < means.shape[1] - 1:
+            a.set_xticklabels('')
+    pb.draw()
+    fig.tight_layout(h_pad=.01) # , rect=(0, 0, 1, .95))
+    return fig