1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24 """The grid search algorithm.
25
26 This file is part of the U{minfx optimisation library<https://sourceforge.net/projects/minfx>}.
27 """
28
29
30 from numpy import array, float64, ones, zeros
31
32
33 from minfx.base_classes import print_iter
34 from minfx.constraint_linear import Constraint_linear
35 from minfx.errors import MinfxError
36
37
38 -def grid(func, args=(), num_incs=None, lower=None, upper=None, incs=None, A=None, b=None, l=None, u=None, c=None, sparseness=None, verbosity=0, print_prefix=""):
39 """The grid search algorithm.
40
41 @param func: The target function. This should take the parameter vector as the first argument and return a single float.
42 @type func: function
43 @keyword args: A tuple of arguments to pass to the function, if needed.
44 @type args: tuple
45 @keyword num_incs: The number of linear increments to be used in the grid search. The length should be equal to the number of parameters and each element corresponds to the number of increments for the respective parameter. This is overridden if the incs argument is supplied.
46 @type num_incs: list of int
47 @keyword lower: The list of lower bounds for the linear grid search. This must be supplied if incs is not.
48 @type lower: list of float
49 @keyword upper: The list of upper bounds for the linear grid search. This must be supplied if incs is not.
50 @type upper: list of float
51 @keyword incs: The parameter increment values. This overrides the num_incs, lower, and upper arguments used in generating a linear grid.
52 @type incs: list of lists
53 @keyword A: The linear constraint matrix A, such that A.x >= b.
54 @type A: numpy rank-2 array
55 @keyword b: The linear constraint scalar vectors, such that A.x >= b.
56 @type b: numpy rank-1 array
57 @keyword l: The lower bound constraint vector, such that l <= x <= u.
58 @type l: list of float
59 @keyword u: The upper bound constraint vector, such that l <= x <= u.
60 @type u: list of float
61 @keyword c: A user supplied constraint function.
62 @type c: function
63 @keyword sparseness: An optional argument for defining sparsity, or regions of the grid to not search over. This is a list whereby each element is a list of two parameters indices. These two parameters will then be assumed to be decoupled and the grid search space between them will be skipped. Symmetry is not observed, so if [0, 1] is sent in, then maybe [1, 0] should be as well.
64 @type sparseness: list of list of int
65 @keyword verbosity: The verbosity level. 0 corresponds to no output, 1 is standard, and higher values cause greater and greater amount of output.
66 @type verbosity: int
67 @keyword print_prefix: The text to place before the printed output.
68 @type print_prefix: str
69 @return: The optimisation information including the parameter vector at the best grid point, the function value at this grid point, the number of iterations (equal to the number of function calls), and a warning.
70 @rtype: tuple of numpy rank-1 array, float, int, str
71 """
72
73
74 if num_incs == None and incs == None:
75 raise MinfxError("Either the incs arg or the num_incs, lower, and upper args must be supplied.")
76 elif num_incs:
77
78 if not isinstance(num_incs, list):
79 raise MinfxError("The num_incs argument '%s' must be a list." % num_incs)
80 if not isinstance(lower, list):
81 raise MinfxError("The lower argument '%s' must be a list." % lower)
82 if not isinstance(upper, list):
83 raise MinfxError("The upper argument '%s' must be a list." % upper)
84
85
86 if len(num_incs) != len(lower):
87 raise MinfxError("The '%s' num_incs and '%s' lower arguments are of different lengths" % (num_incs, lower))
88 if len(num_incs) != len(upper):
89 raise MinfxError("The '%s' num_incs and '%s' upper arguments are of different lengths" % (num_incs, upper))
90
91
92 if num_incs == [] or incs == []:
93
94 if verbosity:
95 print("Cannot run a grid search on a model with zero parameters, directly calculating the function value.")
96
97
98 x0 = zeros(0, float64)
99
100
101 fk = func(x0)
102
103
104 return x0, fk, 1, "No optimisation"
105
106
107
108 if num_incs:
109 n = len(num_incs)
110 else:
111 n = len(incs)
112 grid_size = 0
113 total_steps = 1
114 step_num = ones(n, int)
115 min_step = ones(n, int)
116 params = zeros((n), float64)
117 min_params = zeros((n), float64)
118 sparseness_flag = False
119 if sparseness != None and len(sparseness):
120 sparseness_flag = True
121
122
123
124 if num_incs:
125 incs = []
126
127
128 for k in range(n):
129 params[k] = lower[k]
130 min_params[k] = lower[k]
131 total_steps = total_steps * num_incs[k]
132 incs.append([])
133
134
135 for i in range(num_incs[k]):
136
137 if num_incs[k] == 1:
138 incs[k].append((lower[k] + upper[k]) / 2.0)
139
140
141 else:
142 incs[k].append(lower[k] + i * (upper[k] - lower[k]) / (num_incs[k] - 1))
143
144
145 else:
146 for k in range(n):
147 total_steps = total_steps * len(incs[k])
148
149
150 if verbosity:
151 if verbosity >= 2:
152 print(print_prefix)
153 print(print_prefix)
154 print(print_prefix + "Grid search")
155 print(print_prefix + "~~~~~~~~~~~")
156
157
158 if A != None and b != None:
159 constraint_flag = 1
160 constraint_linear = Constraint_linear(A, b)
161 c = constraint_linear.func
162 if verbosity >= 3:
163 print(print_prefix + "Linear constraint matrices.")
164 print(print_prefix + "A: " + repr(A))
165 print(print_prefix + "b: " + repr(b))
166
167
168 elif l != None and u != None:
169 constraint_flag = 1
170 raise MinfxError("Bound constraints are not implemented yet.")
171
172
173 elif c != None:
174 constraint_flag = 1
175
176
177 else:
178 constraint_flag = 0
179
180
181 f_min = 1e300
182
183
184 if verbosity:
185 print("\n" + print_prefix + "Searching through %s grid nodes." % total_steps)
186
187
188 if total_steps >= 1e8:
189 raise MinfxError("A grid search of size %s is too large." % total_steps)
190
191
192 k = 0
193 curr_dim = 0
194 for i in range(total_steps):
195
196 skip = False
197
198
199 if constraint_flag:
200 ci = c(params)
201 if min(ci) < 0.0:
202 if verbosity >= 3:
203 print_iter(k, min_params, print_prefix=print_prefix)
204 print(print_prefix + "Constraint violated, skipping grid point.")
205 print(print_prefix + "ci: " + repr(ci))
206 print("")
207 skip = True
208
209
210 if sparseness_flag:
211
212 for block in sparseness:
213
214 if curr_dim == block[0]:
215 if step_num[block[1]] != min_step[block[1]]:
216 skip = True
217 break
218
219
220 min_found = False
221 if not skip:
222
223 f = func(*(params,)+args)
224
225
226 if f < f_min:
227 f_min = f
228 min_params = 1.0 * params
229 min_step = 1.0 * step_num
230 min_found = True
231
232
233 if verbosity:
234 print_iter(k=k, xk=min_params, fk=f_min, print_prefix=print_prefix)
235
236
237 grid_size = grid_size + 1
238
239
240 if verbosity >= 2:
241 if not min_found:
242 print_iter(k=k, xk=min_params, fk=f, print_prefix=print_prefix)
243 if verbosity >= 3:
244 if min_found:
245 print(print_prefix + "Minimum found.")
246 print(print_prefix + "%-20s%-20s" % ("Increment:", step_num))
247 print(print_prefix + "%-20s%-20s" % ("Current dimension:", curr_dim))
248 print(print_prefix + "%-20s%-20s" % ("Params:", params))
249 print(print_prefix + "%-20s%-20s" % ("Min params:", min_params))
250 print(print_prefix + "%-20s%-20s" % ("Min indices:", min_step))
251 print(print_prefix + "%-20s%-20.8g" % ("f:", f))
252 print(print_prefix + "%-20s%-20.8g\n" % ("Min f:", f_min))
253
254
255 k = k + 1
256
257
258 for j in range(n):
259 if step_num[j] < len(incs[j]):
260 step_num[j] = step_num[j] + 1
261 params[j] = incs[j][step_num[j]-1]
262 if j > curr_dim:
263 curr_dim = j
264 break
265 else:
266 step_num[j] = 1
267 params[j] = incs[j][0]
268
269
270 return min_params, f_min, grid_size, None
271
272
274 """The grid search algorithm.
275
276 @param func: The target function. This should take the parameter vector as the first argument and return a single float.
277 @type func: function
278 @keyword args: A tuple of arguments to pass to the function, if needed.
279 @type args: tuple
280 @keyword points: The array of grid points to search over.
281 @type points: list or array of lists of float
282 @keyword verbosity: The verbosity level. 0 corresponds to no output, 1 is standard, and higher values cause greater and greater amount of output.
283 @type verbosity: int
284 @keyword print_prefix: The text to place before the printed output.
285 @type print_prefix: str
286 @return: The optimisation information including the parameter vector at the best grid point, the function value at this grid point, the number of iterations (equal to the number of function calls), and a warning.
287 @rtype: tuple of numpy rank-1 array, float, int, str
288 """
289
290
291 total_steps = len(points)
292 n = len(points[0])
293
294
295 f_min = 1e300
296
297
298 if verbosity:
299 print("\n" + print_prefix + "Searching through %s grid nodes." % total_steps)
300
301
302 if total_steps >= 1e8:
303 raise MinfxError("A grid search of size %s is too large." % total_steps)
304
305
306 for k in range(total_steps):
307
308 f = func(*(points[k],)+args)
309
310
311 if f < f_min:
312
313 f_min = f
314 min_params = 1.0 * points[k]
315
316
317 if verbosity:
318 print_iter(k=k, xk=min_params, fk=f_min, print_prefix=print_prefix)
319
320
321 if verbosity >= 2:
322 if f != f_min:
323 print_iter(k=k, xk=min_params, fk=f, print_prefix=print_prefix)
324 if verbosity >= 3:
325 print(print_prefix + "%-20s%-20s" % ("Params:", repr(points[k])))
326 print(print_prefix + "%-20s%-20s" % ("Min params:", repr(min_params)))
327 print(print_prefix + "%-20s%-20g\n" % ("f:", f))
328 print(print_prefix + "%-20s%-20g\n" % ("Min f:", f_min))
329
330
331 return min_params, f_min, total_steps, None
332
333
334 -def grid_split(divisions=None, lower=None, upper=None, inc=None, A=None, b=None, l=None, u=None, c=None):
335 """Generator method yielding arrays of grid points.
336
337 This method will loop over the grid points one-by-one, generating a list of points and yielding these for each subdivision.
338
339
340 @keyword divisions: The number of grid subdivisions.
341 @type divisions: int
342 @keyword lower: The lower bounds of the grid search which must be equal to the number of parameters in the model.
343 @type lower: array of numbers
344 @keyword upper: The upper bounds of the grid search which must be equal to the number of parameters in the model.
345 @type upper: array of numbers
346 @keyword inc: The increments for each dimension of the space for the grid search. The number of elements in the array must equal to the number of parameters in the model.
347 @type inc: array of int
348 @keyword A: The linear constraint matrix A, such that A.x >= b.
349 @type A: numpy rank-2 array
350 @keyword b: The linear constraint scalar vectors, such that A.x >= b.
351 @type b: numpy rank-1 array
352 @keyword l: The lower bound constraint vector, such that l <= x <= u.
353 @type l: list of float
354 @keyword u: The upper bound constraint vector, such that l <= x <= u.
355 @type u: list of float
356 @keyword c: A user supplied constraint function.
357 @type c: function
358 @return: A list of grid points for each subdivision is yielded.
359 @rtype: list of list of float
360 """
361
362
363 n = len(inc)
364
365
366 if A != None and b != None:
367 constraint_flag = True
368 constraint_linear = Constraint_linear(A, b)
369 c = constraint_linear.func
370 if verbosity >= 3:
371 print(print_prefix + "Linear constraint matrices.")
372 print(print_prefix + "A: " + repr(A))
373 print(print_prefix + "b: " + repr(b))
374
375
376 elif l != None and u != None:
377 constraint_flag = True
378 raise MinfxError("Bound constraints are not implemented yet.")
379
380
381 elif c != None:
382 constraint_flag = True
383
384
385 else:
386 constraint_flag = False
387
388
389 total_pts = 1
390 for i in range(n):
391 total_pts = total_pts * inc[i]
392
393
394 indices = zeros(n, int)
395 pts = zeros((total_pts, n), float64)
396
397
398 incs = []
399 for k in range(n):
400 incs.append([])
401 for i in range(inc[k]):
402 incs[k].append(lower[k] + i * (upper[k] - lower[k]) / (inc[k] - 1))
403
404
405 for i in range(total_pts):
406
407 for j in range(n):
408
409 pts[i][j] = incs[j][indices[j]]
410
411
412 for j in range(n):
413 if indices[j] < inc[j]-1:
414 indices[j] += 1
415 break
416 else:
417 indices[j] = 0
418
419
420 if constraint_flag:
421
422 pts_trimmed = []
423 for i in range(total_pts):
424
425 ci = c(pts[i])
426
427
428 if min(ci) >= 0.0:
429 pts_trimmed.append(pts[i])
430
431
432 else:
433 pts_trimmed = pts
434
435
436 pts_trimmed = array(pts_trimmed)
437
438
439 total_pts = len(pts_trimmed)
440
441
442 size_float = total_pts / float(divisions)
443 size = int(size_float)
444 if size_float % 1:
445 size = size + 1
446
447
448 for i in range(min(divisions, total_pts)):
449
450 start = i * size
451
452
453 if i != divisions - 1:
454 end = (i + 1) * size
455 else:
456 end = total_pts
457
458
459 yield pts[start: end]
460