From de201e10072e0fe712a78ba5609cec0b99d7d7fd Mon Sep 17 00:00:00 2001 From: Johann Dreo Date: Thu, 6 Sep 2012 20:52:15 +0200 Subject: [PATCH 1/2] add an evaluator wrapper that keep the best individual found so far --- eo/src/eo | 1 + eo/src/eoEvalKeepBest.h | 107 +++++++++++++++++++++++++++++++++++ eo/test/CMakeLists.txt | 1 + eo/test/t-eoEvalKeepBest.cpp | 78 +++++++++++++++++++++++++ 4 files changed, 187 insertions(+) create mode 100644 eo/src/eoEvalKeepBest.h create mode 100644 eo/test/t-eoEvalKeepBest.cpp diff --git a/eo/src/eo b/eo/src/eo index fd660b7b9..b87ca57f1 100644 --- a/eo/src/eo +++ b/eo/src/eo @@ -78,6 +78,7 @@ #include #include #include +#include // Continuators - all include eoContinue.h #include diff --git a/eo/src/eoEvalKeepBest.h b/eo/src/eoEvalKeepBest.h new file mode 100644 index 000000000..0287cf1de --- /dev/null +++ b/eo/src/eoEvalKeepBest.h @@ -0,0 +1,107 @@ + +/* + This library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; + version 2 of the License. + + This library is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public + License along with this library; if not, write to the Free Software + Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + + © 2012 Thales group + + Authors: + Johann Dreo +*/ + +#ifndef eoEvalKeepBest_H +#define eoEvalKeepBest_H + +#include +#include + +/** + Evaluate with the given evaluator and keep the best individual found so far. + + This is useful if you use a non-monotonic algorithm, such as CMA-ES, where the + population's best fitness can decrease between two generations. This is + sometime necessary and one can't use elitist replacors, as one do not want to + introduce a bias in the population. + + The eoEvalBestKeep is a wrapper around a classical evaluator, that keep the + best individual it has found since its instanciation. + + To get the best individual, you have to call best_element() on the + eoEvalKeepBest itself, and not on the population (or else you would get the + best individual found at the last generation). + + Example: + + MyEval true_eval; + eoEvalKeepBest wrapped_eval( true_eval ); + + // as an interesting side effect, you will get the best individual since + // initalization. + eoPop pop( my_init ); + eoPopLoopEval loop_eval( wrapped_eval ); + + loop_eval( pop ); + + eoEasyEA algo( …, wrapped_eval, … ); + + algo(pop); + + // do not use pop.best_element()! + std::cout << wrapped_eval.best_element() << std::endl; + + @ingroup Evaluation + */ +template class eoEvalKeepBest : public eoEvalFunc, public eoValueParam +{ + public : + eoEvalKeepBest(eoEvalFunc& _func, std::string _name = "VeryBest. ") + : eoValueParam(EOT(), _name), func(_func) {} + + virtual void operator()(EOT& sol) + { + if( sol.invalid() ) { + func(sol); // evaluate + + // if there is no best kept + if( this->value().invalid() ) { + // take the first individual as best + this->value() = sol; + } else { + // if sol is better than the kept individual + if( sol.fitness() > this->value().fitness() ) { + this->value() = sol; + } + } + } // if invalid + } + + //! Return the best individual found so far. + EOT best_element() + { + return this->value(); + } + + /** Reset the best individual to the given one. If no individual is + * provided, the next evaluated one will be taken as a reference. + */ + void reset( const EOT& new_best = EOT() ) + { + this->value() = new_best; + } + + protected : + eoEvalFunc& func; +}; + +#endif diff --git a/eo/test/CMakeLists.txt b/eo/test/CMakeLists.txt index 5f098d3bf..db415bebc 100644 --- a/eo/test/CMakeLists.txt +++ b/eo/test/CMakeLists.txt @@ -34,6 +34,7 @@ ENDIF() ###################################################################################### SET (TEST_LIST + t-eoEvalKeepBest t-eoInitVariableLength t-eofitness t-eoRandom diff --git a/eo/test/t-eoEvalKeepBest.cpp b/eo/test/t-eoEvalKeepBest.cpp new file mode 100644 index 000000000..fdcdae866 --- /dev/null +++ b/eo/test/t-eoEvalKeepBest.cpp @@ -0,0 +1,78 @@ +#include + +#include +#include +#include +#include "real_value.h" + +using namespace std; + +int main(int argc, char* argv[]) +{ + + typedef eoReal EOT; + + eoParser parser(argc, argv); // for user-parameter reading + eoState state; // keeps all things allocated + + + /********************************************* + * problem or representation dependent stuff * + *********************************************/ + + // The evaluation fn - encapsulated into an eval counter for output + eoEvalFuncPtr&> + main_eval( real_value ); // use a function defined in real_value.h + + // wrap the evaluation function in a call counter + eoEvalFuncCounter eval_counter(main_eval); + + // the genotype - through a genotype initializer + eoRealInitBounded& init = make_genotype(parser, state, EOT()); + + // Build the variation operator (any seq/prop construct) + eoGenOp& op = make_op(parser, state, init); + + + /********************************************* + * Now the representation-independent things * + *********************************************/ + + + // initialize the population - and evaluate + // yes, this is representation indepedent once you have an eoInit + eoPop& pop = make_pop(parser, state, init); + + // stopping criteria + eoContinue & term = make_continue(parser, state, eval_counter); + + // things that are called at each generation + eoCheckPoint & checkpoint = make_checkpoint(parser, state, eval_counter, term); + + // wrap the evaluator in another one that will keep the best individual + // evaluated so far + eoEvalKeepBest eval_keep( eval_counter ); + + // algorithm + eoAlgo& ea = make_algo_scalar(parser, state, eval_keep, checkpoint, op); + + + /*************************************** + * Now, call functors and DO something * + ***************************************/ + + // to be called AFTER all parameters have been read! + make_help(parser); + + // evaluate intial population AFTER help and status in case it takes time + apply(eval_keep, pop); + + std::clog << "Best individual after initialization and " << eval_counter.value() << " evaluations" << std::endl; + std::cout << eval_keep.best_element() << std::endl; + + ea(pop); // run the ea + + std::cout << "Best individual after search and " << eval_counter.value() << " evaluations" << std::endl; + // you can also call value(), because it is an eoParam + std::cout << eval_keep.value() << std::endl; +} From 2c2e9d0ca3d75db83f88df9b67ba180bbf854d0b Mon Sep 17 00:00:00 2001 From: Johann Dreo Date: Tue, 18 Sep 2012 16:56:38 +0200 Subject: [PATCH 2/2] better variance computation, use Knuth online robust algorithm, add a test for variance computation --- edo/CMakeLists.txt | 2 +- edo/src/edoEstimatorNormalMono.h | 112 ++++++++++++++++++------------- edo/test/CMakeLists.txt | 1 + edo/test/t-variance.cpp | 38 +++++++++++ 4 files changed, 106 insertions(+), 47 deletions(-) create mode 100644 edo/test/t-variance.cpp diff --git a/edo/CMakeLists.txt b/edo/CMakeLists.txt index 3cb182b3d..9a07cd619 100644 --- a/edo/CMakeLists.txt +++ b/edo/CMakeLists.txt @@ -125,7 +125,7 @@ SET(SAMPLE_SRCS) ADD_SUBDIRECTORY(src) ADD_SUBDIRECTORY(application) -#ADD_SUBDIRECTORY(test) +ADD_SUBDIRECTORY(test) ADD_SUBDIRECTORY(doc) ###################################################################################### diff --git a/edo/src/edoEstimatorNormalMono.h b/edo/src/edoEstimatorNormalMono.h index 12db5188f..9e80145d3 100644 --- a/edo/src/edoEstimatorNormalMono.h +++ b/edo/src/edoEstimatorNormalMono.h @@ -39,64 +39,84 @@ Authors: template < typename EOT > class edoEstimatorNormalMono : public edoEstimator< edoNormalMono< EOT > > { -public: - typedef typename EOT::AtomType AtomType; - - class Variance - { public: - Variance() : _sumvar(0){} + typedef typename EOT::AtomType AtomType; - void update(AtomType v) + //! Knuth's algorithm, online variance, numericably stable + class Variance { - _n++; + public: + Variance() : _n(0), _mean(0), _M2(0) {} - AtomType d = v - _mean; - - _mean += 1 / _n * d; - _sumvar += (_n - 1) / _n * d * d; - } - - AtomType get_mean() const {return _mean;} - AtomType get_var() const {return _sumvar / (_n - 1);} - AtomType get_std() const {return sqrt( get_var() );} - - private: - AtomType _n; - AtomType _mean; - AtomType _sumvar; - }; - -public: - edoNormalMono< EOT > operator()(eoPop& pop) - { - unsigned int popsize = pop.size(); - assert(popsize > 0); - - unsigned int dimsize = pop[0].size(); - assert(dimsize > 0); - - std::vector< Variance > var( dimsize ); - - for (unsigned int i = 0; i < popsize; ++i) - { - for (unsigned int d = 0; d < dimsize; ++d) + void update(AtomType x) { - var[d].update( pop[i][d] ); + _n++; + + AtomType delta = x - _mean; + + _mean += delta / _n; + _M2 += delta * ( x - _mean ); + } + + AtomType mean() const {return _mean;} + + //! Population variance + AtomType var_n() const { + assert( _n > 0 ); + return _M2 / _n; + } + + /** Sample variance (using Bessel's correction) + * is an unbiased estimate of the population variance, + * but it has uniformly higher mean squared error + */ + AtomType var() const { + assert( _n > 1 ); + return _M2 / (_n - 1); + } + + //! Population standard deviation + AtomType std_n() const {return sqrt( var_n() );} + + //! Sample standard deviation, is a biased estimate of the population standard deviation + AtomType std() const {return sqrt( var() );} + + private: + AtomType _n; + AtomType _mean; + AtomType _M2; + }; + + public: + edoNormalMono< EOT > operator()(eoPop& pop) + { + unsigned int popsize = pop.size(); + assert(popsize > 0); + + unsigned int dimsize = pop[0].size(); + assert(dimsize > 0); + + std::vector< Variance > var( dimsize ); + + for (unsigned int i = 0; i < popsize; ++i) + { + for (unsigned int d = 0; d < dimsize; ++d) + { + var[d].update( pop[i][d] ); } } - EOT mean( dimsize ); - EOT variance( dimsize ); + EOT mean( dimsize ); + EOT variance( dimsize ); - for (unsigned int d = 0; d < dimsize; ++d) + for (unsigned int d = 0; d < dimsize; ++d) { - mean[d] = var[d].get_mean(); - variance[d] = var[d].get_var(); + mean[d] = var[d].mean(); + variance[d] = var[d].var_n(); } - return edoNormalMono< EOT >( mean, variance ); - } + return edoNormalMono< EOT >( mean, variance ); + } }; #endif // !_edoEstimatorNormalMono_h diff --git a/edo/test/CMakeLists.txt b/edo/test/CMakeLists.txt index 7a3129a1d..ffc0cfdb8 100644 --- a/edo/test/CMakeLists.txt +++ b/edo/test/CMakeLists.txt @@ -34,6 +34,7 @@ INCLUDE_DIRECTORIES(${CMAKE_SOURCE_DIR}/application/common) SET(SOURCES #t-cholesky + t-variance t-edoEstimatorNormalMulti t-mean-distance t-bounderno diff --git a/edo/test/t-variance.cpp b/edo/test/t-variance.cpp new file mode 100644 index 000000000..ca6c12261 --- /dev/null +++ b/edo/test/t-variance.cpp @@ -0,0 +1,38 @@ + +#include +#include + +#include +#include +#include + +int main() +{ + typedef eoReal Vec; + + eoPop pop; + for( unsigned int i=1; i<7; ++i) { + Vec indiv(1,i); + pop.push_back( indiv ); + std::clog << indiv << " "; + } + std::clog << std::endl; + + edoEstimatorNormalMono estimator; + + edoNormalMono distrib = estimator(pop); + + Vec ex_mean(1,3.5); + Vec ex_var(1,17.5/6); + Vec es_mean = distrib.mean(); + Vec es_var = distrib.variance(); + + std::cout << "expected mean=" << ex_mean << " variance=" << ex_var << std::endl; + std::cout << "estimated mean=" << es_mean << " variance=" << es_var << std::endl; + + for( unsigned int i=0; i