plot_regression_benchmarks.py 10.9 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
#!/usr/bin/env python
# encoding: utf-8
#| Copyright Inria May 2015
#| This project has received funding from the European Research Council (ERC) under
#| the European Union's Horizon 2020 research and innovation programme (grant
#| agreement No 637972) - see http://www.resibots.eu
#|
#| Contributor(s):
#|   - Jean-Baptiste Mouret (jean-baptiste.mouret@inria.fr)
#|   - Antoine Cully (antoinecully@gmail.com)
#|   - Konstantinos Chatzilygeroudis (konstantinos.chatzilygeroudis@inria.fr)
#|   - Federico Allocati (fede.allocati@gmail.com)
#|   - Vaios Papaspyros (b.papaspyros@gmail.com)
#|   - Roberto Rama (bertoski@gmail.com)
#|
#| This software is a computer library whose purpose is to optimize continuous,
#| black-box functions. It mainly implements Gaussian processes and Bayesian
#| optimization.
#| Main repository: http://github.com/resibots/limbo
#| Documentation: http://www.resibots.eu/limbo
#|
#| This software is governed by the CeCILL-C license under French law and
#| abiding by the rules of distribution of free software.  You can  use,
#| modify and/ or redistribute the software under the terms of the CeCILL-C
#| license as circulated by CEA, CNRS and INRIA at the following URL
#| "http://www.cecill.info".
#|
#| As a counterpart to the access to the source code and  rights to copy,
#| modify and redistribute granted by the license, users are provided only
#| with a limited warranty  and the software's author,  the holder of the
#| economic rights,  and the successive licensors  have only  limited
#| liability.
#|
#| In this respect, the user's attention is drawn to the risks associated
#| with loading,  using,  modifying and/or developing or reproducing the
#| software by the user in light of its specific status of free software,
#| that may mean  that it is complicated to manipulate,  and  that  also
#| therefore means  that it is reserved for developers  and  experienced
#| professionals having in-depth computer knowledge. Users are therefore
#| encouraged to load and test the software's suitability as regards their
#| requirements in conditions enabling the security of their systems and/or
#| data to be ensured and,  more generally, to use and operate it in the
#| same conditions as regards security.
#|
#| The fact that you are presently reading this means that you have had
#| knowledge of the CeCILL-C license and that you accept its terms.
47
#|
48
49
from glob import glob
from collections import defaultdict, OrderedDict
50
51
52
53
from datetime import datetime
import platform
import multiprocessing
import os
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96

try:
    from waflib import Logs
    def print_log(c, s): Logs.pprint(c, s)
except: # not in waf
    def print_log(c, s): print(s)

try:
    import numpy as np
    from pylab import *
    import brewer2mpl
    bmap = brewer2mpl.get_map('Set2', 'qualitative', 8)
    colors = bmap.mpl_colors
    plot_ok = True
except:
    plot_ok = False
    Logs.pprint('YELLOW', 'WARNING: numpy/matplotlib not found: no plot of the BO benchmark results')

params = {
    'axes.labelsize' : 12,
    'text.fontsize' : 12,
    'axes.titlesize': 12,
    'legend.fontsize' : 12,
    'xtick.labelsize': 12,
    'ytick.labelsize' : 12,
    'figure.figsize' : [8, 8]
}
rcParams.update(params)

def load_data():
    files = glob("regression_benchmark_results/*/*/*.dat")
    points = defaultdict(lambda : defaultdict(lambda : defaultdict(dict)))
    times_learn = defaultdict(lambda : defaultdict(lambda : defaultdict(dict)))
    times_query = defaultdict(lambda : defaultdict(lambda : defaultdict(dict)))
    mses = defaultdict(lambda : defaultdict(lambda : defaultdict(dict)))
    for f in files:
        if 'test.dat' not in f and 'data.dat' not in f:
            fs = f.split("/")
            func, exp, bench = fs[-1][:-4], fs[-2], fs[-3]
            var = 'limbo'
            if func[-4:] == '_gpy':
                func = func[:-4]
                var = 'GPy'
97
98
            if func[-6:] == '_libgp':
                func = func[:-6]
99
                var = 'libGP'
100
            exp = exp[4:]
101

102
103
            text_file = open(f, "r")
            txt_d = text_file.readlines()
104
105
106
            n_models = int(txt_d[0].strip().split()[-1])
            var_base = var
            for i in range(0, len(txt_d), n_models+1):
107
108
109
                line = txt_d[i].strip().split()
                dim = int(line[0])
                pts = int(line[1])
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130

                for j in range(0, n_models):
                    if n_models > 1:
                        var = var_base + '-' + str(j+1)
                    line2 = txt_d[i+j+1].strip().split()
                    time_learn = float(line2[0])
                    time_query = float(line2[1])
                    mse = float(line2[2])
                    if len(line2) == 4:
                        var = var_base + '-' + line2[-1]
                    # print(bench,var,func,dim)
                    if not (var in points[bench][func][dim]):
                        points[bench][func][dim][var] = []
                        times_learn[bench][func][dim][var] = []
                        times_query[bench][func][dim][var] = []
                        mses[bench][func][dim][var] = []

                    points[bench][func][dim][var].append(pts)
                    times_learn[bench][func][dim][var].append(time_learn)
                    times_query[bench][func][dim][var].append(time_query)
                    mses[bench][func][dim][var].append(mse)
131
132
133
134
135
136
137
138
139
140
    return points,times_learn,times_query,mses

def custom_ax(ax):
    ax.spines['right'].set_visible(False)
    ax.spines['top'].set_visible(False)
    ax.get_xaxis().tick_bottom()
    ax.get_yaxis().tick_left()
    ax.set_axisbelow(True)
    ax.grid(axis='x', color="0.9", linestyle='-')

141
def plot_ax(ax, data, points, labely, disp_legend=True, disp_xaxis=False):
142
    labels = []
143
    replicates = 0
144
    kk = 0
145
146
147
148

    # sort the variants
    sorted_keys = sorted(points.keys(), reverse=True)

149
    # for each variant
150
151
    for vv in range(len(sorted_keys)):
        var = sorted_keys[vv]
152
153
        labels.append(var)
        var_p = points[var]
154
        var_data = data[var]
155
156
157
158
159
160

        pp = {}

        for i in range(len(var_p)):
            if var_p[i] not in pp:
                pp[var_p[i]] = []
161
            pp[var_p[i]].append(var_data[i])
162
163
164
165
166
167
168
169
170
171
        
        pp = OrderedDict(sorted(pp.items()))

        x_axis = pp.keys()
        dd = pp.values()

        y_axis = []
        y_axis_75 = []
        y_axis_25 = []
        for i in range(len(dd)):
172
            replicates = len(dd[i])
173
174
175
176
177
            y_axis.append(np.median(dd[i]))
            y_axis_75.append(np.percentile(dd[i], 75))
            y_axis_25.append(np.percentile(dd[i], 25))

        c_kk = colors[kk%len(colors)]
178
        ax.plot(x_axis, y_axis, '-o', color=c_kk, linewidth=3, markersize=5)
179
180
181
        ax.fill_between(x_axis, y_axis_75, y_axis_25, color=c_kk, alpha=0.15, linewidth=2)
        kk = kk + 1
    
182
183
184
185
    if disp_legend:
        ax.legend(labels)
    if disp_xaxis:
        ax.set_xlabel('Number of points')
186
187
    ax.set_ylabel(labely)
    custom_ax(ax)
188

189
190
191
    return replicates

def plot_data(bench, func, dim, mses, query_times, learning_times, points, rst_file):
192
193
194
    name = func+'_'+str(dim)
    fig, ax = plt.subplots(3, sharex=True)

195
    replicates = plot_ax(ax[0], mses, points, 'Mean Squared Error')
196
197
198
    plot_ax(ax[1], query_times, points, 'Querying time in ms', False)
    plot_ax(ax[2], learning_times, points, 'Learning time in seconds', False, True)

199
    fig.tight_layout()
200
    fig.savefig('regression_benchmark_results/'+bench+'_figs/'+name+'.png')
201
202
    close()

203
204
205
206
207
208
    rst_file.write(func.title() + " in " + str(dim) + "D\n")
    rst_file.write("-----------------\n\n")
    rst_file.write(str(replicates) + " replicates \n\n")
    rst_file.write(".. figure:: fig_benchmarks/" + bench + "_figs/" + name + ".png\n\n")

def plot(points,times_learn,times_query,mses,rst_file):
209
210
    # for each benchmark configuration
    for bench in points.keys():
211
212
213
214
215
216
        fig_dir = os.getcwd() + '/regression_benchmark_results/'+bench+'_figs/'
        try:
            os.makedirs(fig_dir)
            print("created :" + fig_dir)
        except:
            print('WARNING: directory \'%s\' could not be created! (it probably exists already)' % fig_dir)
217
218
219
220
221
        # for each function
        for func in points[bench].keys():
            # for each dimension
            for dim in points[bench][func].keys():
                print('plotting for benchmark: ' + bench + ', the function: ' + func + ' for dimension: ' + str(dim))
222
223
                name = bench+'_'+func+'_'+str(dim)

224
                plot_data(bench, func, dim, mses[bench][func][dim], times_query[bench][func][dim], times_learn[bench][func][dim], points[bench][func][dim], rst_file)
225
226
227
228
229

def plot_all():
    if not plot_ok:
        print_log('YELLOW', "No plot")
        return
230
231
232
233
234
235
236
237

    rst_file = open("regression_benchmark_results/regression_benchmarks.rst", "w")
    rst_file.write("Gaussian process regression benchmarks\n")
    rst_file.write("======================================\n\n")
    date = "{:%B %d, %Y}".format(datetime.datetime.now())
    node = platform.node()
    rst_file.write("*" + date + "* -- " + node + " (" + str(multiprocessing.cpu_count()) + " cores)\n\n")

238
    rst_file.write("- We compare to GPy (https://github.com/SheffieldML/GPy) and libGP (https://github.com/mblum/libgp)\n")
239
240
241
242
    rst_file.write("- Mean Squared Error: lower is better\n")
    rst_file.write("- Learning time: lower is better\n")
    rst_file.write("- Querying time (for 1 by 1 query points): lower is better\n\n")
    rst_file.write("- In each replicate, all the variants (see below for the variants available) are using exactly the same data\n\n")
243
    rst_file.write("- The data are uniformly sampled and some noise is added (according to the variance of the data).\n\n")
244
245
246
247

    rst_file.write("Naming convention\n")
    rst_file.write("------------------\n\n")

248
249
    rst_file.write("- GP-SE-Full-Rprop: Limbo with Squared Exponential kernel where the signal noise, signal variance and kernel lengthscales are optimized via Maximum Likelihood Estimation with the Rprop optimizer (default for limbo)\n")
    rst_file.write("- GP-SE-Rprop: Limbo with Squared Exponential kernel where the signal variance and kernel lengthscales are optimized via Maximum Likelihood Estimation with the Rprop optimizer (default for limbo) and where the signal noise is not optimized but set to a default value: 0.01\n")
250
    rst_file.write("- libGP-SE-Full: libGP with Squared Exponential kernel where the signal noise, signal variance and kernel lengthscales are optimized via Maximum Likelihood Estimation with the Rprop optimizer (the only one that libGP has)\n")
251
    rst_file.write("- GPy: GPy with Squared Exponential kernel where the signal noise, signal variance and kernel lengthscales are optimized via Maximum Likelihood Estimation (with the L-BFGS-B optimizer --- `check scipy.optimize.fmin_l_bfgs_b= <https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.fmin_l_bfgs_b.html>`_)\n\n")
252

253
254
255
    print('loading data...')
    points,times_learn,times_query,mses = load_data()
    print('data loaded')
256
    plot(points,times_learn,times_query,mses,rst_file)
257
258
259

if __name__ == "__main__":
    plot_all()