GridSearch.h
Go to the documentation of this file.
1 //===========================================================================
2 /*!
3  *
4  *
5  * \brief GridSearch.h
6  *
7  *
8  *
9  * \author O. Krause
10  * \date 2010
11  *
12  *
13  * \par Copyright 1995-2015 Shark Development Team
14  *
15  * <BR><HR>
16  * This file is part of Shark.
17  * <http://image.diku.dk/shark/>
18  *
19  * Shark is free software: you can redistribute it and/or modify
20  * it under the terms of the GNU Lesser General Public License as published
21  * by the Free Software Foundation, either version 3 of the License, or
22  * (at your option) any later version.
23  *
24  * Shark is distributed in the hope that it will be useful,
25  * but WITHOUT ANY WARRANTY; without even the implied warranty of
26  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
27  * GNU Lesser General Public License for more details.
28  *
29  * You should have received a copy of the GNU Lesser General Public License
30  * along with Shark. If not, see <http://www.gnu.org/licenses/>.
31  *
32  */
33 //===========================================================================
34 
35 #ifndef SHARK_ALGORITHMS_GRIDSEARCH_H
36 #define SHARK_ALGORITHMS_GRIDSEARCH_H
37 
38 
40 #include <shark/Rng/GlobalRng.h>
41 #include <boost/foreach.hpp>
42 
43 #include <boost/serialization/vector.hpp>
44 
45 
46 namespace shark {
47 
48 //!
49 //! \brief Optimize by trying out a grid of configurations
50 //!
51 //! \par
52 //! The GridSearch class allows for the definition of a grid in
53 //! parameter space. It does a simple one-step optimization over
54 //! the grid by trying out every possible parameter combination.
55 //! Please note that the computation effort grows exponentially
56 //! with the number of parameters.
57 //!
58 //! \par
59 //! If you only want to try a subset of the grid, consider using
60 //! the PointSearch class instead.
61 //! A more sophisticated (less exhaustive) grid search variant is
62 //! available with the NestedGridSearch class.
63 //!
64 class GridSearch : public AbstractSingleObjectiveOptimizer<RealVector >
65 {
66 public:
68  m_configured=false;
69  }
70 
71  /// \brief From INameable: return the class name.
72  std::string name() const
73  {
74  return "GridSearch";
75  }
76 
77  //! uniform initialization for all parameters
78  //! \param params number of model parameters to optimize
79  //! \param min smallest parameter value
80  //! \param max largest parameter value
81  //! \param numSections total number of values in the interval
82  void configure(size_t params, double min, double max, size_t numSections)
83  {
84  SIZE_CHECK(params>=1);
85  RANGE_CHECK(min<=max);
86  SIZE_CHECK(numSections>=1);
87  m_nodeValues.resize(params);
88  for (size_t i = 0; i < numSections; i++)
89  {
90  double section = min + i * (max - min) / (numSections - 1.0);
91  BOOST_FOREACH(std::vector<double>& node,m_nodeValues)
92  node.push_back(section);
93  }
94  m_configured=true;
95  }
96 
97  //! individual definition for every parameter
98  //! \param min smallest value for every parameter
99  //! \param max largest value for every parameter
100  //! \param sections total number of values for every parameter
101  void configure(const std::vector<double>& min, const std::vector<double>& max, const std::vector<size_t>& sections)
102  {
103  size_t params = min.size();
104  SIZE_CHECK(sections.size() == params);
105  SIZE_CHECK(max.size() == params);
106  RANGE_CHECK(min <= max);
107 
108  m_nodeValues.resize(params);
109  for (size_t dimension = 0; dimension < params; dimension++)
110  {
111  size_t numSections = sections[dimension];
112  std::vector<double>& node = m_nodeValues[dimension];
113  node.resize(numSections);
114 
115  if ( numSections == 1 )
116  {
117  node[0] = (( min[dimension] + max[dimension] ) / 2.0);
118  }
119  else for (size_t section = 0; section < numSections; section++)
120  {
121  node[section] = (min[dimension] + section * (max[dimension] - min[dimension]) / (numSections - 1.0));
122  }
123  }
124  m_configured=true;
125  }
126 
127 
128  //! special case for 2D grid, individual definition for every parameter
129  //! \param min1 smallest value for first parameter
130  //! \param max1 largest value for first parameter
131  //! \param sections1 total number of values for first parameter
132  //! \param min2 smallest value for second parameter
133  //! \param max2 largest value for second parameter
134  //! \param sections2 total number of values for second parameter
135  void configure(double min1, double max1, size_t sections1, double min2, double max2, size_t sections2)
136  {
137  RANGE_CHECK(min1<=max1);
138  RANGE_CHECK(min2<=max2);
139  RANGE_CHECK(sections1 > 0);
140  RANGE_CHECK(sections2 > 0);
141 
142  m_nodeValues.resize(2u);
143 
144  if ( sections1 == 1 ) {
145  m_nodeValues[0].push_back(( min1 + max1 ) / 2.0);
146  } else {
147  for (size_t section = 0; section < sections1; section++)
148  m_nodeValues[0].push_back(min1 + section * (max1 - min1) / (sections1 - 1.0));
149  }
150 
151  if ( sections2 == 1 ) {
152  m_nodeValues[1].push_back(( min2 + max2 ) / 2.0);
153  } else {
154  for (size_t section = 0; section < sections2; section++)
155  m_nodeValues[1].push_back(min2 + section * (max2 - min2) / (sections2 - 1.0));
156  }
157  }
158 
159  //! special case for line search
160  //! \param min1 smallest value for first parameter
161  //! \param max1 largest value for first parameter
162  //! \param sections1 total number of values for first parameter
163  void configure(double min1, double max1, size_t sections1)
164  {
165  RANGE_CHECK(min1<=max1);
166  RANGE_CHECK(sections1 > 0);
167 
168  m_nodeValues.resize(1u);
169 
170  if ( sections1 == 1 ) {
171  m_nodeValues[0].push_back(( min1 + max1 ) / 2.0);
172  } else {
173  for (size_t section = 0; section < sections1; section++)
174  m_nodeValues[0].push_back(min1 + section * (max1 - min1) / (sections1 - 1.0));
175  }
176  }
177 
178 
179  //! uniform definition of the values to test for all parameters
180  //! \param params number of model parameters to optimize
181  //! \param values values used for every coordinate
182  void configure(size_t params, const std::vector<double>& values)
183  {
184  SIZE_CHECK(params > 0);
185  SIZE_CHECK(values.size() > 0);
186  m_nodeValues.resize(params);
187  BOOST_FOREACH(std::vector<double>& node,m_nodeValues)
188  node=values;
189  m_configured=true;
190  }
191 
192  //! individual definition for every parameter
193  //! \param values values used. The first dimension is the parameter, the second dimension is the node.
194  void configure(const std::vector<std::vector<double> >& values)
195  {
196  SIZE_CHECK(values.size() > 0);
197  m_nodeValues = values;
198  m_configured=true;
199  }
200 
201  //from ISerializable
202  virtual void read( InArchive & archive )
203  {
204  archive>>m_nodeValues;
205  archive>>m_configured;
206  archive>>m_best.point;
207  archive>>m_best.value;
208  }
209 
210  virtual void write( OutArchive & archive ) const
211  {
212  archive<<m_nodeValues;
213  archive<<m_configured;
214  archive<<m_best.point;
215  archive<<m_best.value;
216  }
217 
218  /*! If Gridsearch wasn't configured before calling this method, it is default constructed
219  * as a net spanning the range [-1,1] in all dimensions with 5 searchpoints (-1,-0.5,0,0.5,1).
220  * so don't forget to scale the parameter-ranges of the objective function!
221  * The startingPoint can actually be anything, only its dimension has to be correct.
222  */
223  virtual void init(ObjectiveFunctionType & objectiveFunction, SearchPointType const& startingPoint) {
224  objectiveFunction.init();
225  checkFeatures(objectiveFunction);
226 
227  if(!m_configured)
228  configure(startingPoint.size(),-1,1,5);
229  SIZE_CHECK(startingPoint.size() == m_nodeValues.size());
230  m_best.point=startingPoint;
231  }
233  /*! Assign linearly progressing grid values to one certain parameter only.
234  * This is especially useful if one parameter needs special treatment
235  * \param index the index of the parameter to which grid values are assigned
236  * \param noOfSections how many grid points should be assigned to that parameter
237  * \param min smallest value for that parameter
238  * \param max largest value for that parameter */
239  void assignLinearRange(size_t index, size_t noOfSections, double min, double max)
240  {
241  SIZE_CHECK( noOfSections >= 1);
242  RANGE_CHECK( min <= max );
243  SIZE_CHECK( index < m_nodeValues.size() );
244 
245  m_nodeValues[index].clear();
246  if ( noOfSections == 1 ) {
247  m_nodeValues[index].push_back(( min+max) / 2.0);
248  }
249  else {
250  m_nodeValues[index].reserve(noOfSections);
251  for (size_t section = 0; section < noOfSections; section++)
252  m_nodeValues[index].push_back(min + section*( max-min ) / ( noOfSections-1.0 ));
253  }
254  }
255 
256  /*! Set exponentially progressing grid values for one certain parameter only.
257  * This is especially useful if one parameter needs special treatment. The
258  * grid points will be filled with values \f$ factor \cdot expbase ^i \f$,
259  * where i does integer steps between min and max.
260  * \param index the index of the parameter that gets new grid values
261  * \param factor the value that the exponential base grid should be multiplied by
262  * \param exp_base the exponential grid will progress on this base (e.g. 2, 10)
263  * \param min the smallest exponent for exp_base
264  * \param max the largest exponent for exp_base */
265  void assignExponentialRange(size_t index, double factor, double exp_base, int min, int max)
266  {
267  SIZE_CHECK( min <= max );
268  RANGE_CHECK( index < m_nodeValues.size() );
269 
270  m_nodeValues[index].clear();
271  m_nodeValues[index].reserve(max-min);
272  for (int section = 0; section <= (max-min); section++)
273  m_nodeValues[index].push_back( factor * std::pow( exp_base, section+min ));
274  }
275 
276  //! Please note that for the grid search optimizer it does
277  //! not make sense to call step more than once, as the
278  //! solution does not improve iteratively.
279  void step(ObjectiveFunctionType const& objectiveFunction) {
280  size_t dimensions = m_nodeValues.size();
281  std::vector<size_t> index(dimensions, 0);
282  m_best.value = 1e100;
283  RealVector point(dimensions);
284 
285  // loop through all grid points
286  while (true)
287  {
288  // define the parameters
289  for (size_t dimension = 0; dimension < dimensions; dimension++)
290  point(dimension) = m_nodeValues[dimension][index[dimension]];
291 
292  // evaluate the model
293  if (objectiveFunction.isFeasible(point))
294  {
295  double error = objectiveFunction.eval(point);
296 
297 #ifdef SHARK_CV_VERBOSE_1
298  std::cout << "." << std::flush;
299 #endif
300 #ifdef SHARK_CV_VERBOSE
301  std::cout << point << "\t" << error << std::endl;
302 #endif
303  if (error < m_best.value)
304  {
305  m_best.value = error;
306  m_best.point = point; // [TG] swap() solution is out, caused ugly memory bug, I changed this back
307  }
308  }
309 
310  // next index
311  size_t dimension = 0;
312  for (; dimension < dimensions; dimension++)
313  {
314  index[dimension]++;
315  if (index[dimension] < m_nodeValues[dimension].size()) break;
316  index[dimension] = 0;
317  }
318  if (dimension == dimensions) break;
319  }
320 #ifdef SHARK_CV_VERBOSE_1
321  std::cout << std::endl;
322 #endif
323  }
324 
325 protected:
326 
327  //! The array columns contain the grid values for the corresponding parameter axis.
328  std::vector<std::vector<double> > m_nodeValues;
329 
331 };
332 
333 
334 //!
335 //! \brief Nested grid search
336 //!
337 //! \par
338 //! The NestedGridSearch class is an iterative optimizer,
339 //! doing one grid search in every iteration. In every
340 //! iteration, it halves the grid extent doubling the
341 //! resolution in every coordinate.
342 //!
343 //! \par
344 //! Although nested grid search is much less exhaustive
345 //! than standard grid search, it still suffers from
346 //! exponential time and memory complexity in the number
347 //! of variables optimized. Therefore, if the number of
348 //! variables is larger than 2 or 3, consider using the
349 //! CMA instead.
350 //!
351 //! \par
352 //! Nested grid search works as follows: The optimizer
353 //! defined a 5x5x...x5 equi-distant grid (depending on
354 //! the search space dimension) on an initially defined
355 //! search cube. During every grid search iteration,
356 //! the error is computed for all grid points.
357 //!Then the grid is moved
358 //! to the best grid point found so far and contracted
359 //! by a factor of two in each dimension. Each call to
360 //! the optimize() function performs one such step.
361 //!
362 //! \par
363 //! Let N denote the number of parameters to optimize.
364 //! To compute the error landscape at the current
365 //! zoom level, the algorithm has to do \f$ 5^N \f$
366 //! error function evaluations in every iteration.
367 //!
368 //! \par
369 //! The grid is always centered around the best
370 //! solution currently known. If this solution is
371 //! located at the boundary, the landscape may exceed
372 //! the parameter range defined m_minimum and m_maximum.
373 //! These invalid landscape values are not used.
374 //!
376 {
377 public:
378  //! Constructor
380  {
381  m_configured=false;
382  }
383 
384  /// \brief From INameable: return the class name.
385  std::string name() const
386  {
387  return "NestedGridSearch";
388  }
389 
390  //!
391  //! \brief Initialization of the nested grid search.
392  //!
393  //! \par
394  //! The min and max arrays define ranges for every parameter to optimize.
395  //! These ranges are strict, that is, the algorithm will not try values
396  //! beyond the range, even if is finds a boundary minimum.
397  //!
398  //! \param min lower end of the parameter range
399  //! \param max upper end of the parameter range
400  void configure(const std::vector<double>& min, const std::vector<double>& max)
401  {
402  size_t dimensions = min.size();
403  SIZE_CHECK(max.size() == dimensions);
404 
405  m_minimum = min;
406  m_maximum = max;
407 
408  m_stepsize.resize(dimensions);
409  m_best.point.resize(dimensions);
410  for (size_t dimension = 0; dimension < dimensions; dimension++)
411  {
412  m_best.point(dimension)=0.5 *(min[dimension] + max[dimension]);
413  m_stepsize[dimension] = 0.25 * (max[dimension] - min[dimension]);
414  }
415  m_configured=true;
416  }
417 
418  //!
419  //! \brief Initialization of the nested grid search.
420  //!
421  //! \par
422  //! The min and max values define ranges for every parameter to optimize.
423  //! These ranges are strict, that is, the algorithm will not try values
424  //! beyond the range, even if is finds a boundary minimum.
425  //!
426  //! \param parameters number of parameters to optimize
427  //! \param min lower end of the parameter range
428  //! \param max upper end of the parameter range
429  void configure(size_t parameters, double min, double max)
430  {
431  SIZE_CHECK(parameters > 0);
432 
433  m_minimum=std::vector<double>(parameters,min);
434  m_maximum=std::vector<double>(parameters,max);
435  m_stepsize=std::vector<double>(parameters,0.25 * (max - min));
436 
437  m_best.point.resize(parameters);
438 
439  double start=0.5 *(min + max);
440  BOOST_FOREACH(double& value,m_best.point)
441  value=start;
442  m_configured=true;
443  }
444 
445  //from ISerializable
446  virtual void read( InArchive & archive )
447  {
448  archive>>m_minimum;
449  archive>>m_maximum;
450  archive>>m_stepsize;
451  archive>>m_configured;
452  archive>>m_best.point;
453  archive>>m_best.value;
454  }
455 
456  virtual void write( OutArchive & archive ) const
457  {
458  archive<<m_minimum;
459  archive<<m_maximum;
460  archive<<m_stepsize;
461  archive<<m_configured;
462  archive<<m_best.point;
463  archive<<m_best.value;
464  }
465 
466  //! if NestedGridSearch was not configured before this call, it is default initialized ti the range[-1,1] for every parameter
467  virtual void init(ObjectiveFunctionType & objectiveFunction, SearchPointType const& startingPoint) {
468  objectiveFunction.init();
469  checkFeatures(objectiveFunction);
470 
471  if(!m_configured)
472  configure(startingPoint.size(),-1,1);
473  SIZE_CHECK(m_stepsize.size()==startingPoint.size());
474 
475  }
477 
478 
479  //! Every call of the optimization member computes the
480  //! error landscape on the current grid. It picks the
481  //! best error value and zooms into the error landscape
482  //! by a factor of 2.
483  void step(ObjectiveFunctionType const& objectiveFunction) {
484  size_t dimensions = m_stepsize.size();
485  SIZE_CHECK(m_minimum.size() == dimensions);
486  SIZE_CHECK(m_maximum.size() == dimensions);
487 
488  m_best.value = 1e99;
489  std::vector<size_t> index(dimensions,0);
490 
491  RealVector point=m_best.point;
492 
493  // loop through the grid
494  while (true)
495  {
496  // compute the grid point,
497 
498  // set the parameters
499  bool compute=true;
500  for (size_t d = 0; d < dimensions; d++)
501  {
502  point(d) += (index[d] - 2.0) * m_stepsize[d];
503  if (point(d) < m_minimum[d] || point(d) > m_maximum[d])
504  {
505  compute = false;
506  break;
507  }
508  }
509 
510  // evaluate the grid point
511  if (compute && objectiveFunction.isFeasible(point))
512  {
513  double error = objectiveFunction.eval(point);
514 
515  // remember the best solution
516  if (error < m_best.value)
517  {
518  m_best.value = error;
519  m_best.point=point;
520  }
521  }
522  // move to the next grid point
523  size_t d = 0;
524  for (; d < dimensions; d++)
525  {
526  index[d]++;
527  if (index[d] <= 4) break;
528  index[d] = 0;
529  }
530  if (d == dimensions) break;
531  }
532  // decrease the step sizes
533  BOOST_FOREACH(double& step,m_stepsize)
534  step *= 0.5;
535  }
536 
537 protected:
538  //! minimum parameter value to check
539  std::vector<double> m_minimum;
540 
541  //! maximum parameter value to check
542  std::vector<double> m_maximum;
543 
544  //! current step size for every parameter
545  std::vector<double> m_stepsize;
546 
548 };
549 
550 
551 //!
552 //! \brief Optimize by trying out predefined configurations
553 //!
554 //! \par
555 //! The PointSearch class is similair to the GridSearch class
556 //! by the property that it optimizes a model in a single pass
557 //! just trying out a predefined number of parameter configurations.
558 //! The main difference is that every parameter configuration has
559 //! to be explicitly defined. It is not possible to define a set
560 //! of values for every axis; see GridSearch for this purpose.
561 //! Thus, the PointSearch class allows for more flexibility.
562 //!
563 //! If no configure method is called, this class just samples random points.
564 //! They are uniformly distributed in [-1,1].
565 //! parameters^2 points but minimum 20 are sampled in this case.
566 //!
568 {
569 public:
570  //! Constructor
572  m_configured=false;
573  }
574 
575  /// \brief From INameable: return the class name.
576  std::string name() const
577  {
578  return "PointSearch";
579  }
580 
581  //! Initialization of the search points.
582  //! \param points array of points to evaluate
583  void configure(const std::vector<RealVector>& points) {
584  m_points=points;
585  m_configured=true;
586  }
587 
588  //!samples random points in the range [min,max]^parameters
589  void configure(size_t parameters,size_t samples, double min,double max) {
590  RANGE_CHECK(min<=max);
591  m_points.resize(samples);
592  for(size_t sample=0; sample!=samples; ++sample)
593  {
594  m_points[sample].resize(parameters);
595  for(size_t param=0; param!=parameters; ++param)
596  {
597  m_points[sample](param)=Rng::uni(min,max);
598  }
599  }
600  m_configured=true;
601  }
602 
603  virtual void read( InArchive & archive )
604  {
605  archive>>m_points;
606  archive>>m_configured;
607  archive>>m_best.point;
608  archive>>m_best.value;
609  }
610 
611  virtual void write( OutArchive & archive ) const
612  {
613  archive<<m_points;
614  archive<<m_configured;
615  archive<<m_best.point;
616  archive<<m_best.value;
617  }
618 
619  //! If the class wasn't configured before, this method samples random uniform distributed points in [-1,1]^n.
620  void init(ObjectiveFunctionType& objectiveFunction, SearchPointType const& startingPoint) {
621  objectiveFunction.init();
622  checkFeatures(objectiveFunction);
623 
624  if(!m_configured)
625  {
626  size_t parameters=startingPoint.size();
627  size_t samples=std::min(sqr(parameters),(size_t)20);
628  configure(parameters,samples,-1,1);
629  }
630  }
632  //! Please note that for the point search optimizer it does
633  //! not make sense to call step more than once, as the
634  //! solution does not improve iteratively.
635  void step(ObjectiveFunctionType const& objectiveFunction) {
636  size_t numPoints = m_points.size();
637  m_best.value = 1e100;
638  size_t bestIndex=0;
639 
640  // loop through all points
641  for (size_t point = 0; point < numPoints; point++)
642  {
643  // evaluate the model
644  if (objectiveFunction.isFeasible(m_points[point]))
645  {
646  double error = objectiveFunction.eval(m_points[point]);
647  if (error < m_best.value)
648  {
649  m_best.value = error;
650  bestIndex=point;
651  }
652  }
653  }
654  m_best.point=m_points[bestIndex];
655  }
656 
657 protected:
658  //! The array holds one parameter configuration in every column.
659  std::vector<RealVector> m_points;
660 
661  //! verbosity level
663 };
664 
665 
666 }
667 #endif