SHOGUN
v3.2.0
|
00001 /* 00002 * This program is free software; you can redistribute it and/or modify 00003 * it under the terms of the GNU General Public License as published by 00004 * the Free Software Foundation; either version 3 of the License, or 00005 * (at your option) any later version. 00006 * 00007 * Written (W) 2013 Roman Votyakov 00008 * Copyright (C) 2012 Jacob Walker 00009 */ 00010 00011 #include <shogun/modelselection/GradientModelSelection.h> 00012 00013 #ifdef HAVE_NLOPT 00014 00015 #include <shogun/evaluation/GradientResult.h> 00016 #include <shogun/modelselection/ParameterCombination.h> 00017 #include <shogun/modelselection/ModelSelectionParameters.h> 00018 #include <shogun/machine/Machine.h> 00019 #include <nlopt.h> 00020 00021 using namespace shogun; 00022 00023 #ifndef DOXYGEN_SHOULD_SKIP_THIS 00024 00026 struct nlopt_params 00027 { 00029 CMachineEvaluation* machine_eval; 00030 00032 CParameterCombination* current_combination; 00033 00035 CMap<TParameter*, CSGObject*>* parameter_dictionary; 00036 00038 bool print_state; 00039 }; 00040 00051 double nlopt_function(unsigned n, const double* x, double* grad, void* func_data) 00052 { 00053 nlopt_params* params=(nlopt_params*)func_data; 00054 00055 CMachineEvaluation* machine_eval=params->machine_eval; 00056 CParameterCombination* current_combination=params->current_combination; 00057 CMap<TParameter*, CSGObject*>* parameter_dictionary=params->parameter_dictionary; 00058 bool print_state=params->print_state; 00059 00060 index_t offset=0; 00061 00062 // set parameters from vector x 00063 for (index_t i=0; i<parameter_dictionary->get_num_elements(); i++) 00064 { 00065 CMapNode<TParameter*, CSGObject*>* node=parameter_dictionary->get_node_ptr(i); 00066 00067 TParameter* param=node->key; 00068 CSGObject* parent=node->data; 00069 00070 if (param->m_datatype.m_ctype==CT_VECTOR || 00071 param->m_datatype.m_ctype==CT_SGVECTOR) 00072 { 00073 REQUIRE(param->m_datatype.m_length_y, "Parameter vector %s has no " 00074 "length\n", param->m_name) 00075 00076 for (index_t j=0; j<*(param->m_datatype.m_length_y); j++) 00077 { 00078 00079 bool result=current_combination->set_parameter(param->m_name, 00080 (float64_t)x[offset++], parent, j); 00081 REQUIRE(result, "Parameter %s not found in combination tree\n", 00082 param->m_name) 00083 } 00084 } 00085 else 00086 { 00087 bool result=current_combination->set_parameter(param->m_name, 00088 (float64_t)x[offset++], parent); 00089 REQUIRE(result, "Parameter %s not found in combination tree\n", 00090 param->m_name) 00091 } 00092 } 00093 00094 // apply current combination to the machine 00095 CMachine* machine=machine_eval->get_machine(); 00096 current_combination->apply_to_machine(machine); 00097 SG_UNREF(machine); 00098 00099 // evaluate the machine 00100 CEvaluationResult* evaluation_result=machine_eval->evaluate(); 00101 CGradientResult* gradient_result=CGradientResult::obtain_from_generic( 00102 evaluation_result); 00103 SG_UNREF(evaluation_result); 00104 00105 if (print_state) 00106 { 00107 gradient_result->print_result(); 00108 } 00109 00110 // get value of the function, gradients and parameter dictionary 00111 SGVector<float64_t> value=gradient_result->get_value(); 00112 CMap<TParameter*, SGVector<float64_t> >* gradient=gradient_result->get_gradient(); 00113 CMap<TParameter*, CSGObject*>* gradient_dictionary= 00114 gradient_result->get_paramter_dictionary(); 00115 SG_UNREF(gradient_result); 00116 00117 offset=0; 00118 00119 // set derivative for each parameter from parameter dictionary 00120 for (index_t i=0; i<parameter_dictionary->get_num_elements(); i++) 00121 { 00122 CMapNode<TParameter*, CSGObject*>* node=parameter_dictionary->get_node_ptr(i); 00123 00124 SGVector<float64_t> derivative; 00125 00126 for (index_t j=0; j<gradient_dictionary->get_num_elements(); j++) 00127 { 00128 CMapNode<TParameter*, CSGObject*>* gradient_node= 00129 gradient_dictionary->get_node_ptr(j); 00130 00131 if (gradient_node->data==node->data && 00132 !strcmp(gradient_node->key->m_name, node->key->m_name)) 00133 { 00134 derivative=gradient->get_element(gradient_node->key); 00135 } 00136 } 00137 00138 REQUIRE(derivative.vlen, "Can't find gradient wrt %s parameter!\n", 00139 node->key->m_name); 00140 00141 memcpy(grad+offset, derivative.vector, sizeof(double)*derivative.vlen); 00142 00143 offset+=derivative.vlen; 00144 } 00145 00146 SG_UNREF(gradient); 00147 SG_UNREF(gradient_dictionary); 00148 00149 return (double)(SGVector<float64_t>::sum(value)); 00150 } 00151 00152 #endif /* DOXYGEN_SHOULD_SKIP_THIS */ 00153 00154 CGradientModelSelection::CGradientModelSelection() : CModelSelection() 00155 { 00156 init(); 00157 } 00158 00159 CGradientModelSelection::CGradientModelSelection(CMachineEvaluation* machine_eval, 00160 CModelSelectionParameters* model_parameters) 00161 : CModelSelection(machine_eval, model_parameters) 00162 { 00163 init(); 00164 } 00165 00166 CGradientModelSelection::~CGradientModelSelection() 00167 { 00168 } 00169 00170 void CGradientModelSelection::init() 00171 { 00172 m_max_evaluations=1000; 00173 m_grad_tolerance=1e-6; 00174 00175 SG_ADD(&m_grad_tolerance, "gradient_tolerance", "Gradient tolerance", 00176 MS_NOT_AVAILABLE); 00177 SG_ADD(&m_max_evaluations, "max_evaluations", "Maximum number of evaluations", 00178 MS_NOT_AVAILABLE); 00179 } 00180 00181 CParameterCombination* CGradientModelSelection::select_model(bool print_state) 00182 { 00183 if (!m_model_parameters) 00184 { 00185 CMachine* machine=m_machine_eval->get_machine(); 00186 00187 CParameterCombination* current_combination=new CParameterCombination(machine); 00188 SG_REF(current_combination); 00189 00190 if (print_state) 00191 { 00192 SG_PRINT("Initial combination:\n"); 00193 current_combination->print_tree(); 00194 } 00195 00196 // get total length of variables 00197 index_t total_variables=current_combination->get_parameters_length(); 00198 00199 // build parameter->value map 00200 CMap<TParameter*, SGVector<float64_t> >* argument= 00201 new CMap<TParameter*, SGVector<float64_t> >(); 00202 current_combination->build_parameter_values_map(argument); 00203 00204 // unroll current parameter combination into vector 00205 SGVector<double> x(total_variables); 00206 index_t offset=0; 00207 00208 for (index_t i=0; i<argument->get_num_elements(); i++) 00209 { 00210 CMapNode<TParameter*, SGVector<float64_t> >* node=argument->get_node_ptr(i); 00211 memcpy(x.vector+offset, node->data.vector, sizeof(double)*node->data.vlen); 00212 offset+=node->data.vlen; 00213 } 00214 00215 SG_UNREF(argument); 00216 00217 // create nlopt object and choose MMA (Method of Moving Asymptotes) 00218 // optimization algorithm 00219 nlopt_opt opt=nlopt_create(NLOPT_LD_MMA, total_variables); 00220 00221 // create lower bound vector (lb=-inf) 00222 SGVector<double> lower_bound(total_variables); 00223 lower_bound.set_const(1e-6); 00224 00225 // create upper bound vector (ub=inf) 00226 SGVector<double> upper_bound(total_variables); 00227 upper_bound.set_const(HUGE_VAL); 00228 00229 // set upper and lower bound 00230 nlopt_set_lower_bounds(opt, lower_bound.vector); 00231 nlopt_set_upper_bounds(opt, upper_bound.vector); 00232 00233 // set maximum number of evaluations 00234 nlopt_set_maxeval(opt, m_max_evaluations); 00235 00236 // set absolute argument tolearance 00237 nlopt_set_xtol_abs1(opt, m_grad_tolerance); 00238 nlopt_set_ftol_abs(opt, m_grad_tolerance); 00239 00240 // build parameter->sgobject map from current parameter combination 00241 CMap<TParameter*, CSGObject*>* parameter_dictionary= 00242 new CMap<TParameter*, CSGObject*>(); 00243 current_combination->build_parameter_parent_map(parameter_dictionary); 00244 00245 // nlopt parameters 00246 nlopt_params params; 00247 00248 params.current_combination=current_combination; 00249 params.machine_eval=m_machine_eval; 00250 params.print_state=print_state; 00251 params.parameter_dictionary=parameter_dictionary; 00252 00253 // choose evaluation direction (minimize or maximize objective function) 00254 if (m_machine_eval->get_evaluation_direction()==ED_MINIMIZE) 00255 { 00256 if (print_state) 00257 SG_PRINT("Minimizing objective function:\n"); 00258 00259 nlopt_set_min_objective(opt, nlopt_function, ¶ms); 00260 } 00261 else 00262 { 00263 if (print_state) 00264 SG_PRINT("Maximizing objective function:\n"); 00265 00266 nlopt_set_max_objective(opt, nlopt_function, ¶ms); 00267 } 00268 00269 // the minimum objective value, upon return 00270 double minf; 00271 00272 // optimize our function 00273 nlopt_result result=nlopt_optimize(opt, x.vector, &minf); 00274 00275 REQUIRE(result>0, "NLopt failed while optimizing objective function!\n"); 00276 00277 if (print_state) 00278 { 00279 SG_PRINT("Best combination:\n"); 00280 current_combination->print_tree(); 00281 } 00282 00283 // clean up 00284 nlopt_destroy(opt); 00285 SG_UNREF(machine); 00286 SG_UNREF(parameter_dictionary); 00287 00288 return current_combination; 00289 } 00290 else 00291 { 00292 SG_NOTIMPLEMENTED 00293 return NULL; 00294 } 00295 } 00296 00297 #endif /* HAVE_NLOPT */