Commits

Ruben Martinez-Cantin committed ad4546b

Adding missing files. Exceptions in criteria.

Comments (0)

Files changed (7)

cmake/FindPython.cmake

+# This code sets the following variables:
+# PYTHON_EXEC         - path to python executable
+# PYTHON_LIBRARIES    - path to the python library
+# PYTHON_INCLUDE_DIRS - path to where Python.h is found
+# PTYHON_SITE_MODULES - path to site-packages
+# PYTHON_ARCH         - name of architecture to be used for platform-specific
+#                       binary modules
+# PYTHON_VERSION      - version of python
+#
+# You can optionally include the version number when using this package
+# like so:
+#   find_package(python 2.6)
+#
+# This code defines a helper function find_python_module(). It can be used
+# like so:
+#   find_python_module(numpy)
+# If numpy is found, the variable PY_NUMPY contains the location of the numpy
+# module. If the module is required add the keyword "REQUIRED":
+#   find_python_module(numpy REQUIRED)
+
+include(FindPackageHandleStandardArgs)
+
+# allow specifying which Python installation to use
+if (NOT PYTHON_EXEC)
+    set(PYTHON_EXEC $ENV{PYTHON_EXEC})
+endif (NOT PYTHON_EXEC)
+
+if (NOT PYTHON_EXEC)
+    find_program(PYTHON_EXEC "python${Python_FIND_VERSION}"
+        PATHS
+        [HKEY_LOCAL_MACHINE\\SOFTWARE\\Python\\PythonCore\\3.1\\InstallPath]
+        [HKEY_LOCAL_MACHINE\\SOFTWARE\\Python\\PythonCore\\3.0\\InstallPath]
+        [HKEY_LOCAL_MACHINE\\SOFTWARE\\Python\\PythonCore\\2.7\\InstallPath]
+        [HKEY_LOCAL_MACHINE\\SOFTWARE\\Python\\PythonCore\\2.6\\InstallPath]
+        [HKEY_LOCAL_MACHINE\\SOFTWARE\\Python\\PythonCore\\2.5\\InstallPath]
+        DOC "Location of python executable to use")
+endif(NOT PYTHON_EXEC)
+
+# if Python is still not found, return
+if (NOT PYTHON_EXEC)
+    # dummy function
+    function(find_python_module module)
+        return()
+    endfunction(find_python_module)
+    return()
+endif()
+
+# On OS X the python executable might be symlinked to the "real" location
+# of the python executable. The header files and libraries are found relative
+# to that path.
+# For CMake 2.6 and below, the REALPATH option is included in the ABSOLUTE option
+if (${CMAKE_MAJOR_VERSION}.${CMAKE_MINOR_VERSION} GREATER 2.6)
+  get_filename_component(PYTHON_EXEC_ "${PYTHON_EXEC}" REALPATH)
+else()
+  get_filename_component(PYTHON_EXEC_ "${PYTHON_EXEC}" ABSOLUTE)
+endif()
+set(PYTHON_EXEC "${PYTHON_EXEC_}" CACHE FILEPATH "Path to Python interpreter")
+
+string(REGEX REPLACE "/bin/python.*" "" PYTHON_PREFIX "${PYTHON_EXEC_}")
+string(REGEX REPLACE "/bin/python.*" "" PYTHON_PREFIX2 "${PYTHON_EXEC}")
+
+execute_process(COMMAND "${PYTHON_EXEC}" "-c"
+    "import sys; print('%d.%d' % (sys.version_info[0],sys.version_info[1]))"
+    OUTPUT_VARIABLE PYTHON_VERSION
+    OUTPUT_STRIP_TRAILING_WHITESPACE)
+string(REPLACE "." "" PYTHON_VERSION_NO_DOTS ${PYTHON_VERSION})
+
+find_library(PYTHON_LIBRARIES 
+    NAMES "python${PYTHON_VERSION_NO_DOTS}" "python${PYTHON_VERSION}"
+    PATHS
+        "${PYTHON_PREFIX}/lib"
+        [HKEY_LOCAL_MACHINE\\SOFTWARE\\Python\\PythonCore\\${PYTHON_VERSION}\\InstallPath]/libs
+    PATH_SUFFIXES "" "python${PYTHON_VERSION}/config"
+    DOC "Python libraries" NO_DEFAULT_PATH)
+
+find_path(PYTHON_INCLUDE_DIRS "Python.h"
+    PATHS
+        "${PYTHON_PREFIX}/include"
+        [HKEY_LOCAL_MACHINE\\SOFTWARE\\Python\\PythonCore\\${PYTHON_VERSION}\\InstallPath]/include
+    PATH_SUFFIXES python${PYTHON_VERSION} python${PYTHON_VERSION}m
+    DOC "Python include directories" NO_DEFAULT_PATH)
+
+execute_process(COMMAND "${PYTHON_EXEC}" "-c"
+    "from distutils.sysconfig import get_python_lib; print(get_python_lib())"
+    OUTPUT_VARIABLE PYTHON_SITE_MODULES_
+    OUTPUT_STRIP_TRAILING_WHITESPACE)
+string(REGEX REPLACE "^${PYTHON_PREFIX2}/" "" PYTHON_SITE_MODULES "${PYTHON_SITE_MODULES_}")
+
+function(find_python_module module)
+    string(TOUPPER ${module} module_upper)
+    if(ARGC GREATER 1)
+        if (ARGV1 STREQUAL "REQUIRED")
+            set(PY_${module}_FIND_REQUIRED TRUE)
+        else()
+            if (ARGV1 STREQUAL "QUIET")
+                set(PY_${module}_FIND_QUIETLY TRUE)
+            endif()
+        endif()
+    endif()
+    if(NOT PY_${module_upper})
+        # A module's location is usually a directory, but for binary modules
+        # it's a .so file.
+        execute_process(COMMAND "${PYTHON_EXEC}" "-c"
+            "import re, ${module}; print(re.compile('/__init__.py.*').sub('',${module}.__file__))"
+            RESULT_VARIABLE _${module}_status
+            OUTPUT_VARIABLE _${module}_location
+            ERROR_QUIET OUTPUT_STRIP_TRAILING_WHITESPACE)
+        if(NOT _${module}_status)
+            set(PY_${module_upper} ${_${module}_location} CACHE STRING
+                "Location of Python module ${module}")
+        endif(NOT _${module}_status)
+    endif(NOT PY_${module_upper})
+    find_package_handle_standard_args(PY_${module} DEFAULT_MSG PY_${module_upper})
+endfunction(find_python_module)
+
+set(PYTHON_ARCH "unknown")
+if(APPLE)
+    set(PYTHON_ARCH "darwin")
+else(APPLE)
+    if(CMAKE_SYSTEM_NAME STREQUAL "Linux")
+        if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64")
+            set(PYTHON_ARCH "linux2")
+        else(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64")
+            set(PYTHON_ARCH "linux")
+        endif(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64")
+    else(CMAKE_SYSTEM_NAME STREQUAL "Linux")
+        if(CMAKE_SYSTEM_NAME STREQUAL "Windows")
+            set(PYTHON_ARCH "windows")
+        endif(CMAKE_SYSTEM_NAME STREQUAL "Windows")
+    endif(CMAKE_SYSTEM_NAME STREQUAL "Linux")
+endif(APPLE)
+
+find_package_handle_standard_args(Python DEFAULT_MSG
+    PYTHON_LIBRARIES PYTHON_INCLUDE_DIRS PYTHON_SITE_MODULES)

include/conditionalbayesprocess.hpp

+/** \file conditionalbayesprocess.hpp
+    \brief Implementes a empirical Bayesian nonparametric process with a 
+    ML, MAP or similar estimate of kernel parameters. */
+/*
+-------------------------------------------------------------------------
+   This file is part of BayesOpt, an efficient C++ library for 
+   Bayesian optimization.
+
+   Copyright (C) 2011-2013 Ruben Martinez-Cantin <rmcantin@unizar.es>
+ 
+   BayesOpt is free software: you can redistribute it and/or modify it 
+   under the terms of the GNU General Public License as published by
+   the Free Software Foundation, either version 3 of the License, or
+   (at your option) any later version.
+
+   BayesOpt 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 General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with BayesOpt.  If not, see <http://www.gnu.org/licenses/>.
+------------------------------------------------------------------------
+*/
+
+
+#ifndef  _CONDITIONAL_BAYES_PROCESS_HPP_
+#define  _CONDITIONAL_BAYES_PROCESS_HPP_
+
+#include "kernelregressor.hpp"
+#include "inneroptimization.hpp"
+
+namespace bayesopt
+{
+
+  /** \addtogroup  NonParametricProcesses */
+  /**@{*/
+
+
+  /**
+   * \brief Empirical Bayesian NonParametric process.
+   */
+  class ConditionalBayesProcess: public KernelRegressor, RBOptimizable
+  {
+  public:
+    ConditionalBayesProcess(size_t dim, bopt_params parameters,
+			  const Dataset& data, randEngine& eng);
+    virtual ~ConditionalBayesProcess();
+
+    /** 
+     * \brief Function that returns the prediction of the GP for a query point
+     * in the hypercube [0,1].
+     * 
+     * @param query in the hypercube [0,1] to evaluate the Gaussian process
+     * @return pointer to the probability distribution.
+     */	
+    virtual ProbabilityDistribution* prediction(const vectord &query) = 0;
+		 		 
+    /** 
+     * \brief Updates the kernel parameters acording with a point
+     * estimate (ML, MAP, etc.)
+     */
+    void updateKernelParameters();
+
+    /** 
+     * \brief Computes the score (eg:likelihood) of the kernel
+     * parameters.  
+     * Warning: To evaluate the score, it is necessary to change the parameters
+     * @param x set of parameters.  
+     * @return score
+     */
+    double evaluate(const vectord &x);
+
+    /** 
+     * \brief Computes the score (eg:likelihood) of the current kernel
+     * parameters.
+     * @param query set of parameters.
+     * @return score
+     */
+    double evaluateKernelParams();
+
+
+  protected:
+    /** 
+     * \brief Computes the negative log likelihood of the data for all
+     * the parameters.
+     * @return value negative log likelihood
+     */
+    virtual double negativeTotalLogLikelihood() = 0;
+
+
+    /** 
+     * \brief Computes the negative log likelihood of the data for the
+     * kernel hyperparameters.
+     * @return value negative log likelihood
+     */
+    virtual double negativeLogLikelihood() = 0;
+
+  private:
+    /**
+     * Computes the negative score of the data using cross validation.
+     * @return negative score
+     */
+    double negativeCrossValidation();
+
+  private:
+    NLOPT_Optimization* kOptimizer;
+  };
+
+
+
+  inline double ConditionalBayesProcess::evaluate(const vectord& x)
+  { 
+    mKernel.setHyperParameters(x);
+    return evaluateKernelParams();
+  };
+
+
+  /**@}*/
+  
+} //namespace bayesopt
+
+#endif
+

include/criteria_atomic.hpp

       mExp = 1;
     };
 
-    int setParameters(const vectord &params)
-    {
-      mExp = static_cast<size_t>(params(0));
-      return 0;
-    };
+    void setParameters(const vectord &params)
+    { mExp = static_cast<size_t>(params(0)); };
 
     size_t nParameters() {return 1;};
 
       mBias = 0.01;
       mExp = 1;
     };
-    int setParameters(const vectord &params)
+
+    void setParameters(const vectord &params)
     {
       mExp = static_cast<size_t>(params(0));
       mBias = params(1);
-      return 0;
     };
 
     size_t nParameters() {return 2;};
       mProc = proc;
       mBeta = 1.0;
     };
-    int setParameters(const vectord &params)
-    {
-      mBeta = params(0);
-      return 0;
-    };
+    void setParameters(const vectord &params)
+    { mBeta = params(0); };
 
     size_t nParameters() {return 1;};
 
       mProc = proc;
       mEpsilon = 0.01;
     };
-    int setParameters(const vectord &params)
-    {
-      mEpsilon = params(0);
-      return 0;
-    };
+    void setParameters(const vectord &params)
+    { mEpsilon = params(0); };
 
     size_t nParameters() {return 1;};
 
   {
   public:
     virtual ~GreedyAOptimality(){};
-    int setParameters(const vectord &params) { return 0; };
+    void setParameters(const vectord &params) {};
     size_t nParameters() {return 0;};
     double operator()( const vectord &x)
     { return -mProc->prediction(x)->getStd(); };
   {
   public:
     virtual ~ExpectedReturn(){};
-    int setParameters(const vectord &params) { return 0; };
+    void setParameters(const vectord &params) { };
     size_t nParameters() {return 0;};
     double operator()( const vectord &x)
     { return mProc->prediction(x)->getMean(); };
   public:
     OptimisticSampling() {};
     virtual ~OptimisticSampling(){};
-    int setParameters(const vectord &params) { return 0; };
+    void setParameters(const vectord &params) {};
     size_t nParameters() {return 0;};
     double operator()( const vectord &x)
     {
   public:
     ThompsonSampling() {};
     virtual ~ThompsonSampling(){};
-    int setParameters(const vectord &params) { return 0; };
+    void setParameters(const vectord &params) { };
     size_t nParameters() {return 0;};
     double operator()( const vectord &x)
     {
       reset();
     };
 
-    int setParameters(const vectord &params)
-    {
-      mExp = static_cast<size_t>(params(0));
-      return 0;
-    };
+    void setParameters(const vectord &params)
+    { mExp = static_cast<size_t>(params(0)); };
 
     size_t nParameters() {return 1;};
     void reset() { nCalls = 0; mExp = 10;};
       reset();
     };
 
-    int setParameters(const vectord &params)
-    {
-      mCoef = params(0);
-      return 0;
-    };
+    void setParameters(const vectord &params)
+    { mCoef = params(0); };
+
     size_t nParameters() {return 1;};
     void reset() { nCalls = 0; mCoef = 5.0;};
     double operator()( const vectord &x)
     unsigned int nCalls;
   };
 
+
   /**
    * \brief Distance in input space. Can be combined with other
    * critera to trade off large changes in input space.
       mW = 1;
     };
     virtual ~InputDistance(){};
-    int setParameters(const vectord &params)
-    {
-      mW = params(0);
-      return 0;
-    };
+    void setParameters(const vectord &params)
+    { mW = params(0); };
     size_t nParameters() {return 1;};
  
     double operator()(const vectord &x)

include/criteria_combined.hpp

       mCriteriaList = list;
     };
 
-    int setParameters(const vectord &theta) 
+    void setParameters(const vectord &theta) 
     {
       using boost::numeric::ublas::subrange;
       const size_t np = mCriteriaList.size();
       if (theta.size() != norm_1(sizes))
 	{
 	  FILE_LOG(logERROR) << "Wrong number of criteria parameters"; 
-	  return -1; 
+	  throw std::invalid_argument("Wrong number of criteria parameters");
 	}
 
       size_t start = 0;
 	  mCriteriaList[i]->setParameters(subrange(theta,start,start+sizes(i)));
 	  start += sizes(i);
 	}
-      return 0;
     };
 
     size_t nParameters() 

include/criteria_functors.hpp

 
     double evaluate(const vectord &x) {return (*this)(x);}
     virtual double operator()(const vectord &x) = 0;
+
     virtual std::string name() = 0;
-    virtual int setParameters(const vectord &params) = 0;
+    virtual void setParameters(const vectord &params) = 0;
     virtual size_t nParameters() = 0;
 
-    //Dummy functions.
+    //Dummy functions. Not all criteria support these methods.
     virtual void reset() { assert(false); };
     void setRandomEngine(randEngine& eng){ mtRandom = &eng; }
 
   };
 
 
-  template <typename CriteriaType> Criteria * create_func()
-  {
-    return new CriteriaType();
-  }
-
-
   /** 
    * \brief Factory model for criterion functions
    * This factory is based on the libgp library by Manuel Blum

src/conditionalbayesprocess.cpp

+
+/*
+-------------------------------------------------------------------------
+   This file is part of BayesOpt, an efficient C++ library for 
+   Bayesian optimization.
+
+   Copyright (C) 2011-2013 Ruben Martinez-Cantin <rmcantin@unizar.es>
+ 
+   BayesOpt is free software: you can redistribute it and/or modify it 
+   under the terms of the GNU General Public License as published by
+   the Free Software Foundation, either version 3 of the License, or
+   (at your option) any later version.
+
+   BayesOpt 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 General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with BayesOpt.  If not, see <http://www.gnu.org/licenses/>.
+------------------------------------------------------------------------
+*/
+
+
+#include "conditionalbayesprocess.hpp"
+#include "log.hpp"
+//#include "optimizekernel.hpp"	
+
+
+namespace bayesopt
+{
+  ConditionalBayesProcess::ConditionalBayesProcess(size_t dim, bopt_params parameters, 
+						   const Dataset& data, randEngine& eng):
+    KernelRegressor(dim,parameters,data,eng)
+  { 
+    // if (mLearnType == L_BAYES)
+    //   {
+    // 	FILE_LOG(logERROR) << "Empirical Bayes model and full Bayes learning are incompatible.";
+    // 	throw std::invalid_argument("Trying full Bayes learning for an empirical Bayes model.");
+    //   }
+
+    size_t nhp = mKernel.nHyperParameters();
+    kOptimizer = new NLOPT_Optimization(this,nhp);
+
+    //TODO: Generalize
+    if (parameters.sc_type == SC_ML)
+      {
+	kOptimizer->setAlgorithm(BOBYQA);    // local search to avoid underfitting
+      }
+    else
+      {
+	kOptimizer->setAlgorithm(COMBINED);
+      }
+    // kOptimizer->setLimits(svectord(nhp,1e-10),svectord(nhp,100.));
+    kOptimizer->setLimits(svectord(nhp,-6),svectord(nhp,6));
+  }
+
+  ConditionalBayesProcess::~ConditionalBayesProcess()
+  {
+    delete kOptimizer;
+  }
+
+
+  void ConditionalBayesProcess::updateKernelParameters()
+  {
+    // if (mLearnType == L_FIXED)
+    //   {
+    // 	FILE_LOG(logDEBUG) << "Fixed hyperparameters. Not learning";
+    //   }
+    // else
+    //   {
+	vectord optimalTheta = mKernel.getHyperParameters();
+	
+	FILE_LOG(logDEBUG) << "Initial kernel parameters: " << optimalTheta;
+	kOptimizer->run(optimalTheta);
+	mKernel.setHyperParameters(optimalTheta);
+	FILE_LOG(logDEBUG) << "Final kernel parameters: " << optimalTheta;	
+      // }
+  };
+
+  double ConditionalBayesProcess::evaluateKernelParams()
+  { 
+    switch(mScoreType)
+      {
+      case SC_MTL:
+	return negativeTotalLogLikelihood();
+      case SC_ML:
+	return negativeLogLikelihood();
+      case SC_MAP:
+	// It is a minus because the prior is the positive and we want
+	// the negative.
+	return negativeLogLikelihood()-mKernel.kernelLogPrior();
+      case SC_LOOCV:
+	return negativeCrossValidation(); 
+      default:
+	FILE_LOG(logERROR) << "Learning type not supported";
+	throw std::invalid_argument("Learning type not supported");
+      }	  
+  }
+
+
+  double ConditionalBayesProcess::negativeCrossValidation()
+  {
+    // This is highly ineffient implementation for comparison purposes.
+    Dataset data(mData);
+
+    size_t n = data.getNSamples();
+    size_t last = n-1;
+    int error = 0;
+    double sum = 0.0;
+
+    matrixd tempF(mMean.mFeatM);
+
+
+    // We take the first element, use it for validation and then paste
+    // it at the end. Thus, after every iteration, the first element
+    // is different and, at the end, all the elements should have
+    // rotated.
+    for(size_t i = 0; i<n; ++i)
+      {
+	// Take the first element
+	const double y = data.getSampleY(0);
+	const vectord x = data.getSampleX(0);
+
+	// Remove it for cross validation
+	data.mX.erase(data.mX.begin()); 
+	utils::erase(data.mY,data.mY.begin());
+	utils::erase_column(mMean.mFeatM,0);
+
+	// Compute the cross validation
+	computeCholeskyCorrelation();
+	precomputePrediction(); 
+	ProbabilityDistribution* pd = prediction(x);
+	sum += log(pd->pdf(y));
+
+	//Paste it back at the end
+	data.addSample(x,y);
+	mMean.mFeatM.resize(mMean.mFeatM.size1(),mMean.mFeatM.size2()+1);  
+	mMean.mFeatM = tempF;
+      }
+    std::cout << "End" << data.getNSamples();
+    return -sum;   //Because we are minimizing.
+  }
+
+} // namespace bayesopt

src/criteria_functors.cpp

 
 namespace bayesopt
 {
+
+  template <typename CriteriaType> Criteria * create_func()
+  {
+    return new CriteriaType();
+  }
+
+
   CriteriaFactory::CriteriaFactory()
   {
     registry["cEI"] = & create_func<ExpectedImprovement>;
     registry["cHedgeRandom"] = & create_func<GP_Hedge_Random>;
   }
 
-  /// Factory method for criterion functions.
-  // Criteria* CriteriaFactory::create(criterium_name name,
-  // 				    KernelRegressor* proc)
-  // {
-  //   Criteria* c_ptr;
-  //   std::vector<Criteria*> list;
-  //   switch (name)
-  //     {
-  //     case C_EI:     c_ptr = new ExpectedImprovement(); break;
-  //     case C_EI_A:   c_ptr = new AnnealedExpectedImprovement(); break;
-  //     case C_LCB:    c_ptr = new LowerConfidenceBound(); break;
-  //     case C_LCB_A:  c_ptr = new AnnealedLowerConfindenceBound(); break;
-  //     case C_POI:    c_ptr = new ProbabilityOfImprovement(); break;
-  //     case C_GREEDY_A_OPTIMALITY: c_ptr = new GreedyAOptimality(); break;
-  //     case C_EXPECTED_RETURN: c_ptr = new ExpectedReturn(); break;
-  //     case C_OPTIMISTIC_SAMPLING: c_ptr = new OptimisticSampling(); break;
-  //     case C_GP_HEDGE: c_ptr = new GP_Hedge(); break;
-  //     case C_GP_HEDGE_RANDOM: c_ptr = new GP_Hedge_Random(); break;
-  //     default:
-  // 	FILE_LOG(logERROR) << "Error in criterium";
-  // 	return NULL;
-  //     }
-  //   if ((name = C_GP_HEDGE) || (name = C_GP_HEDGE_RANDOM))
-  //     {
-  // 	for(size_t i = 0; i < N_ALGORITHMS_IN_GP_HEDGE; ++i)
-  // 	  {
-  // 	    list.push_back(create(ALGORITHMS_IN_GP_HEDGE[i],proc)); 
-  // 	  }
-  // 	c_ptr->init(proc,list);
-  //     }
-  //   else
-  //     {
-  // 	c_ptr->init(proc);
-  //     }
-  //   return c_ptr;
-  // };
-
 
   /** 
    * \brief Factory model for criterion functions
 	FILE_LOG(logERROR) << "Error: Fatal error while parsing "
 			   << "kernel function: " << os 
 			   << " not found";
-	throw std::invalid_argument("Criteria not found " + os);
+	throw std::invalid_argument("Parsing error: Criteria not found: " + os);
 	return NULL;
       } 
     cFunc = it->second();