In [1]:
import itertools
import os.path
from sets import Set
import sys
sys.path.append('../plume/')

import tables
from matplotlib.gridspec import GridSpec, GridSpecFromSubplotSpec, SubplotSpec
from neatplots.tools import HSplitAxes, SharedAxesGrid
from neatplots.predefined import four as palette

from plume.config import instantiate

%pylab inline

import latexstyle
latexstyle.setup()

styles = [dict(marker=m, c=c) for m, c in zip(['o', '^', 'v', 's'], palette.thin[1:] + [palette.thin[0]])]
gray = (0.3, 0.3, 0.3)


Populating the interactive namespace from numpy and matplotlib

In [2]:
error_type_labels = {
    'rmse': 'RMISE',
    'wrmse': 'WRIMSE',
    'ise': 'ISE',
    'wise': 'WISE',
    'reward': 'QRSim reward',
    'log_likelihood': 'log likelihood'
}
kern_labels = {
    'RBFKernel': r'$k_{\mathrm{SE}}(r)$',
    'Matern32Kernel': r'$k_{5/2}(r)$',
    'Matern52Kernel': r'$k_{3/2}(r)$',
    'ExponentialKernel': r'$k_{\mathrm{exp}}(r)$'
}
scenario_labels = {
    'SingleSource_noiseless': 'G-NF-SS',
    'SingleSource_dispersion_noiseless': 'D-NF-SS',
    'MultipleSource_dispersion_noiseless': 'D-NF-MS',
    'SingleSource_dispersion': 'D-SN-SS',
    'MultipleSource_dispersion': 'D-SN-MS'
}

In [3]:
def get_filename(kern, scenario):
    return '../Data/Kernels/{kern}__{scenario}_'.format(kern=kern, scenario=scenario)

def load_scenario_data(
        scenario, kernels=['RBFKernel', 'Matern52Kernel', 'Matern32Kernel', 'ExponentialKernel']):
    return [tables.open_file(get_filename(kern, scenario)) for kern in kernels
            if os.path.isfile(get_filename(kern, scenario))]

def plot_error_vs_lengthscale(ax, data, error_type, **kwargs):
    kern_name = kern_labels[data.root.conf[0]['kernel'][1]]
    
    zp_err = data.get_node('/' + error_type + '/zero_pred', 'value').read()
    err = data.get_node('/' + error_type + '/gp_pred', 'value').read()
    percentage = err / zp_err
    mean = np.mean(percentage, axis=2)
    std_err = np.std(percentage, axis=2) / np.sqrt(percentage.shape[2])
    print mean.min()

    ax.errorbar(
        data.root.lengthscales.read(), mean[:, 0], std_err[:, 0], fmt='o-', label=kern_name, **kwargs)
    
def plot_error_vs_variance(ax, data, error_type, **kwargs):
    kern_name = kern_labels[data.root.conf[0]['kernel'][1]]
    
    zp_err = data.get_node('/' + error_type + '/zero_pred', 'value').read()
    err = data.get_node('/' + error_type + '/gp_pred', 'value').read()
    percentage = err / zp_err
    mean = np.mean(percentage, axis=2)
    std_err = np.std(percentage, axis=2) / np.sqrt(percentage.shape[2])
    ax.errorbar(
        data.root.variances.read(), mean[0, :], std_err[0, :], fmt='o-', label=kern_name, **kwargs)
    
def create_kernel_plot_matrix(scenarios, error_types, plot_fn=plot_error_vs_lengthscale):
    fig = plt.figure(figsize=(5.35, 6))
    outlier_plot_ratio = 0.4

    class KernelErrorPlot(SharedAxesGrid):
        def _create_axes(self, subplot_spec, sharex, sharey):
            return HSplitAxes(fig, subplot_spec, outlier_plot_ratio, sharex=sharex, sharey=sharey)
        
        def _get_share_xy(self, split):
            return split.bottom, split.both
        
        def _plot(self, split, error_type, scenario):
            for kern, style in zip(load_scenario_data(scenario), styles):
                plot_fn(split.both, kern, error_type, **style)
            split.both.axhline(1.0, c=gray)
    
    legend_height = 0.05
    grid = GridSpec(
        2, 1, height_ratios=(1.0 - legend_height, legend_height), hspace=2 * legend_height,
        bottom=0.05, left=0.1, right=0.975, top=0.95)
            
    p = KernelErrorPlot(error_types, scenarios, SubplotSpec(grid, 0), wspace=0.15, hspace=0.15)
    p.axes.top.loglog()
    p.axes.bottom.semilogx()
    
    latexstyle.style_axes(p.axes.both)
    p.axes.top.xaxis.set_ticks_position('none')
    plt.setp(p.axes.top.get_xticklabels(), visible=False)
    for ax in p.axes_by_row[:-1]:
        plt.setp(ax.bottom.get_xticklabels(), visible=False)
    for ax in p.axes_by_col[1:]:
        plt.setp(ax.both.get_yticklabels(), visible=False)
    
    p.axes.bottom.set_ylim([0, 1.5])
    p.axes.top.set_ylim(bottom=1.85)
    
    #p.axes.top.spines.get('bottom').set_visible(False)
    
    p.axes_by_row[-1].set_xlabel(r'$\ell$', labelpad=0)
    for ax, error_type in zip(p.axes_by_row, error_types):
        ax.set_ylabel('Normalized ' + error_type_labels[error_type], labelpad=0)
    for ax, scenario in zip(p.axes_by_col, scenarios):
        ax[0].set_title(scenario_labels[scenario])
        
    ax_legend = fig.add_subplot(grid[1])
    ax_legend.set_axis_off()
    ax_legend.legend(
        *p.axes[0].top.get_legend_handles_labels(), ncol=4, loc='upper center',
        bbox_to_anchor=(0.5, 1.0), frameon=False, columnspacing=1.5, handletextpad=0.2)
    return fig

def create_kernel_plot_matrix_no_outlier(scenarios, error_types, plot_fn=plot_error_vs_lengthscale):
    fig = plt.figure(figsize=(4, 4))
    outlier_plot_ratio = 0.4

    class KernelErrorPlotNoOutlier(SharedAxesGrid):
        def _create_axes(self, subplot_spec, sharex, sharey):
            return fig.add_subplot(subplot_spec, sharex=sharex, sharey=sharey)
        
        def _plot(self, ax, error_type, scenario):
            for kern, style in zip(load_scenario_data(scenario), styles):
                plot_fn(ax, kern, error_type, **style)
    
    legend_height = 0.05
    grid = GridSpec(
        1, 1,
        bottom=0.05, left=0.1, right=0.975, top=0.95)
            
    p = KernelErrorPlotNoOutlier(error_types, scenarios, SubplotSpec(grid, 0), wspace=0.15, hspace=0.15)
    p.axes.semilogx()
    
    latexstyle.style_axes(p.axes)
    for ax in p.axes_by_row[:-1]:
        plt.setp(ax.get_xticklabels(), visible=False)
    for ax in p.axes_by_col[1:]:
        plt.setp(ax.get_yticklabels(), visible=False)
    
    p.axes.set_ylim([0, 1.0])
    p.axes.set_xlim([1e-12, 1e-3])

    p.axes_by_row[-1].set_xlabel(r'$\ell$', labelpad=0)
    for ax, error_type in zip(p.axes_by_col[0], error_types):
        ax.set_ylabel('Normalized ' + error_type_labels[error_type], labelpad=0)
    for ax, scenario in zip(p.axes_by_col, scenarios):
        ax[0].set_title(scenario_labels[scenario])
    return fig

In [11]:
fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
data = load_scenario_data('MultipleSource_dispersion_noiseless', ['Matern32Kernel'])[0]
plot_error_vs_lengthscale(ax, data, 'rmse')


0.738989782914

In [12]:
scenarios = [
    'SingleSource_noiseless', 'SingleSource_dispersion_noiseless',
    'MultipleSource_dispersion_noiseless']
error_types = ['rmse', 'wrmse', 'reward']
fig = create_kernel_plot_matrix(scenarios, error_types)


0.0653250415133
0.0241897595197
0.026769946721
0.0463211777615
0.689801095199
0.620667473894
0.614745583718
0.609315586172
0.774182672555
0.747706610512
0.738989782914
0.721272447498
0.0102780555595
0.00353327877041
0.00650658817541
0.0229845960885
0.629322436405
0.600192140331
0.595949456894
0.626102617988
0.697259632999
0.665629237837
0.670629483246
0.719368179655
0.00187003849052
0.000397464518286
0.000628261827373
0.00263469346868
0.269083057158
0.223657752708
0.204336508021
0.186055313331
0.588362681203
0.484093870201
0.459066682783
0.413060720893

In [13]:
fig.savefig('../../thesis/plots/lengthscales.pdf')

In [14]:
scenarios = ['SingleSource_dispersion', 'MultipleSource_dispersion']
error_types = ['rmse', 'wrmse', 'reward']
fig = create_kernel_plot_matrix_no_outlier(scenarios, error_types, plot_error_vs_variance)



In [15]:
fig.savefig('../../thesis/plots/kvar.pdf')

In [4]:
scenarios = [
    'SingleSource_noiseless', 'SingleSource_dispersion_noiseless',
    'MultipleSource_dispersion_noiseless']
error_types = ['log_likelihood']

fig = plt.figure(figsize=(5.35, 1.6))

class KernelErrorPlot(SharedAxesGrid):
    def _create_axes(self, subplot_spec, sharex, sharey):
        return fig.add_subplot(subplot_spec, sharex=sharex, sharey=sharey)
    
    def _plot(self, ax, error_type, scenario):
        for kern, style in zip(load_scenario_data(scenario), styles):
            kern_name = kern_labels[kern.root.conf[0]['kernel'][1]]
    
            ll = kern.get_node('/' + error_type + '/gp_pred', 'value').read()
            mean = np.mean(ll, axis=2)
            std_err = np.std(ll, axis=2) / np.sqrt(ll.shape[2])
            print scenario, kern_name, mean.max(), kern.root.lengthscales.read()[2], mean[2]

            ax.errorbar(
                kern.root.lengthscales.read(), mean[:, 0], std_err[:, 0],
                fmt='o-', label=kern_name, **style)

legend_height = 0.1
grid = GridSpec(
    2, 1, height_ratios=(1.0 - legend_height, legend_height), hspace=0.5,
    bottom=0.1, left=0.1, right=0.975, top=0.875)
        
p = KernelErrorPlot(error_types, scenarios, SubplotSpec(grid, 0))
latexstyle.style_axes(p.axes)

p.axes.set_xlabel(r'$\ell$', labelpad=0)
p.axes.semilogx()
p.axes.set_ylim(-1500, 2600)
p.axes_by_col[0].set_ylabel('log likelihood', labelpad=0)
for ax in p.axes_by_col[1:]:
    plt.setp(ax.get_yticklabels(), visible=False)
for ax, scenario in zip(p.axes_by_col, scenarios):
    ax[0].set_title(scenario_labels[scenario])
    
ax_legend = fig.add_subplot(grid[1])
ax_legend.set_axis_off()
ax_legend.legend(
    *p.axes[0].get_legend_handles_labels(), ncol=4, loc='upper center',
    bbox_to_anchor=(0.5, 1.0), frameon=False, columnspacing=1.5, handletextpad=0.2)


SingleSource_noiseless $k_{\mathrm{SE}}(r)$ -40.8717977747 5.0 [-40.87179777]
SingleSource_noiseless $k_{3/2}(r)$ 1572.14411251 5.0 [-506.11348652]
SingleSource_noiseless $k_{5/2}(r)$ 1663.50187359 5.0 [-605.51904275]
SingleSource_noiseless $k_{\mathrm{exp}}(r)$ 378.563311236 5.0 [-768.99507388]
SingleSource_dispersion_noiseless $k_{\mathrm{SE}}(r)$ -89.0232045161 5.0 [-89.02320452]
SingleSource_dispersion_noiseless $k_{3/2}(r)$ 1639.56902224 5.0 [-501.81205398]
SingleSource_dispersion_noiseless $k_{5/2}(r)$ 2472.04584731 5.0 [-601.83970125]
SingleSource_dispersion_noiseless $k_{\mathrm{exp}}(r)$ 378.684956855 5.0 [-766.41494249]
MultipleSource_dispersion_noiseless $k_{\mathrm{SE}}(r)$ -321.013680009 5.0 [-321.01368001]
MultipleSource_dispersion_noiseless $k_{3/2}(r)$ 1859.67376873 5.0 [-600.47838597]
MultipleSource_dispersion_noiseless $k_{5/2}(r)$ 2359.98191122 5.0 [-671.72938996]
MultipleSource_dispersion_noiseless $k_{\mathrm{exp}}(r)$ 341.14847592 5.0 [-795.53243089]
Out[4]:
<matplotlib.legend.Legend at 0x10854b290>

In [5]:
fig.savefig('../../thesis/plots/loglikelihood.pdf')

In [ ]: