SHOGUN
v2.0.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) 2007-2008 Vojtech Franc 00008 * Written (W) 2007-2009 Soeren Sonnenburg 00009 * Copyright (C) 2007-2009 Fraunhofer Institute FIRST and Max-Planck-Society 00010 */ 00011 00012 #include <shogun/labels/Labels.h> 00013 #include <shogun/mathematics/Math.h> 00014 #include <shogun/lib/Time.h> 00015 #include <shogun/base/Parameter.h> 00016 #include <shogun/base/Parallel.h> 00017 #include <shogun/machine/LinearMachine.h> 00018 #include <shogun/classifier/svm/SVMOcas.h> 00019 #include <shogun/features/DotFeatures.h> 00020 #include <shogun/labels/Labels.h> 00021 #include <shogun/labels/BinaryLabels.h> 00022 00023 using namespace shogun; 00024 00025 CSVMOcas::CSVMOcas() 00026 : CLinearMachine() 00027 { 00028 init(); 00029 } 00030 00031 CSVMOcas::CSVMOcas(E_SVM_TYPE type) 00032 : CLinearMachine() 00033 { 00034 init(); 00035 method=type; 00036 } 00037 00038 CSVMOcas::CSVMOcas( 00039 float64_t C, CDotFeatures* traindat, CLabels* trainlab) 00040 : CLinearMachine() 00041 { 00042 init(); 00043 C1=C; 00044 C2=C; 00045 00046 set_features(traindat); 00047 set_labels(trainlab); 00048 } 00049 00050 00051 CSVMOcas::~CSVMOcas() 00052 { 00053 } 00054 00055 bool CSVMOcas::train_machine(CFeatures* data) 00056 { 00057 SG_INFO("C=%f, epsilon=%f, bufsize=%d\n", get_C1(), get_epsilon(), bufsize); 00058 SG_DEBUG("use_bias = %i\n", get_bias_enabled()) ; 00059 00060 ASSERT(m_labels); 00061 ASSERT(m_labels->get_label_type() == LT_BINARY); 00062 if (data) 00063 { 00064 if (!data->has_property(FP_DOT)) 00065 SG_ERROR("Specified features are not of type CDotFeatures\n"); 00066 set_features((CDotFeatures*) data); 00067 } 00068 ASSERT(features); 00069 00070 int32_t num_vec=features->get_num_vectors(); 00071 lab = SGVector<float64_t>(num_vec); 00072 for (int32_t i=0; i<num_vec; i++) 00073 lab[i] = ((CBinaryLabels*)m_labels)->get_label(i); 00074 00075 w=SGVector<float64_t>(features->get_dim_feature_space()); 00076 w.zero(); 00077 00078 if (num_vec!=lab.vlen || num_vec<=0) 00079 SG_ERROR("num_vec=%d num_train_labels=%d\n", num_vec, lab.vlen); 00080 00081 SG_FREE(old_w); 00082 old_w=SG_CALLOC(float64_t, w.vlen); 00083 bias=0; 00084 old_bias=0; 00085 00086 tmp_a_buf=SG_CALLOC(float64_t, w.vlen); 00087 cp_value=SG_CALLOC(float64_t*, bufsize); 00088 cp_index=SG_CALLOC(uint32_t*, bufsize); 00089 cp_nz_dims=SG_CALLOC(uint32_t, bufsize); 00090 cp_bias=SG_CALLOC(float64_t, bufsize); 00091 00092 float64_t TolAbs=0; 00093 float64_t QPBound=0; 00094 int32_t Method=0; 00095 if (method == SVM_OCAS) 00096 Method = 1; 00097 ocas_return_value_T result = svm_ocas_solver( get_C1(), num_vec, get_epsilon(), 00098 TolAbs, QPBound, get_max_train_time(), bufsize, Method, 00099 &CSVMOcas::compute_W, 00100 &CSVMOcas::update_W, 00101 &CSVMOcas::add_new_cut, 00102 &CSVMOcas::compute_output, 00103 &CSVMOcas::sort, 00104 &CSVMOcas::print, 00105 this); 00106 00107 SG_INFO("Ocas Converged after %d iterations\n" 00108 "==================================\n" 00109 "timing statistics:\n" 00110 "output_time: %f s\n" 00111 "sort_time: %f s\n" 00112 "add_time: %f s\n" 00113 "w_time: %f s\n" 00114 "solver_time %f s\n" 00115 "ocas_time %f s\n\n", result.nIter, result.output_time, result.sort_time, 00116 result.add_time, result.w_time, result.qp_solver_time, result.ocas_time); 00117 00118 SG_FREE(tmp_a_buf); 00119 00120 primal_objective = result.Q_P; 00121 00122 uint32_t num_cut_planes = result.nCutPlanes; 00123 00124 SG_DEBUG("num_cut_planes=%d\n", num_cut_planes); 00125 for (uint32_t i=0; i<num_cut_planes; i++) 00126 { 00127 SG_DEBUG("cp_value[%d]=%p\n", i, cp_value); 00128 SG_FREE(cp_value[i]); 00129 SG_DEBUG("cp_index[%d]=%p\n", i, cp_index); 00130 SG_FREE(cp_index[i]); 00131 } 00132 00133 SG_FREE(cp_value); 00134 cp_value=NULL; 00135 SG_FREE(cp_index); 00136 cp_index=NULL; 00137 SG_FREE(cp_nz_dims); 00138 cp_nz_dims=NULL; 00139 SG_FREE(cp_bias); 00140 cp_bias=NULL; 00141 00142 SG_FREE(old_w); 00143 old_w=NULL; 00144 00145 return true; 00146 } 00147 00148 /*---------------------------------------------------------------------------------- 00149 sq_norm_W = sparse_update_W( t ) does the following: 00150 00151 W = oldW*(1-t) + t*W; 00152 sq_norm_W = W'*W; 00153 00154 ---------------------------------------------------------------------------------*/ 00155 float64_t CSVMOcas::update_W( float64_t t, void* ptr ) 00156 { 00157 float64_t sq_norm_W = 0; 00158 CSVMOcas* o = (CSVMOcas*) ptr; 00159 uint32_t nDim = (uint32_t) o->w.vlen; 00160 float64_t* W=o->w.vector; 00161 float64_t* oldW=o->old_w; 00162 00163 for(uint32_t j=0; j <nDim; j++) 00164 { 00165 W[j] = oldW[j]*(1-t) + t*W[j]; 00166 sq_norm_W += W[j]*W[j]; 00167 } 00168 o->bias=o->old_bias*(1-t) + t*o->bias; 00169 sq_norm_W += CMath::sq(o->bias); 00170 00171 return( sq_norm_W ); 00172 } 00173 00174 /*---------------------------------------------------------------------------------- 00175 sparse_add_new_cut( new_col_H, new_cut, cut_length, nSel ) does the following: 00176 00177 new_a = sum(data_X(:,find(new_cut ~=0 )),2); 00178 new_col_H = [sparse_A(:,1:nSel)'*new_a ; new_a'*new_a]; 00179 sparse_A(:,nSel+1) = new_a; 00180 00181 ---------------------------------------------------------------------------------*/ 00182 int CSVMOcas::add_new_cut( 00183 float64_t *new_col_H, uint32_t *new_cut, uint32_t cut_length, 00184 uint32_t nSel, void* ptr) 00185 { 00186 CSVMOcas* o = (CSVMOcas*) ptr; 00187 CDotFeatures* f = o->features; 00188 uint32_t nDim=(uint32_t) o->w.vlen; 00189 float64_t* y = o->lab.vector; 00190 00191 float64_t** c_val = o->cp_value; 00192 uint32_t** c_idx = o->cp_index; 00193 uint32_t* c_nzd = o->cp_nz_dims; 00194 float64_t* c_bias = o->cp_bias; 00195 00196 float64_t sq_norm_a; 00197 uint32_t i, j, nz_dims; 00198 00199 /* temporary vector */ 00200 float64_t* new_a = o->tmp_a_buf; 00201 memset(new_a, 0, sizeof(float64_t)*nDim); 00202 00203 for(i=0; i < cut_length; i++) 00204 { 00205 f->add_to_dense_vec(y[new_cut[i]], new_cut[i], new_a, nDim); 00206 00207 if (o->use_bias) 00208 c_bias[nSel]+=y[new_cut[i]]; 00209 } 00210 00211 /* compute new_a'*new_a and count number of non-zerou dimensions */ 00212 nz_dims = 0; 00213 sq_norm_a = CMath::sq(c_bias[nSel]); 00214 for(j=0; j < nDim; j++ ) { 00215 if(new_a[j] != 0) { 00216 nz_dims++; 00217 sq_norm_a += new_a[j]*new_a[j]; 00218 } 00219 } 00220 00221 /* sparsify new_a and insert it to the last column of sparse_A */ 00222 c_nzd[nSel] = nz_dims; 00223 c_idx[nSel]=NULL; 00224 c_val[nSel]=NULL; 00225 00226 if(nz_dims > 0) 00227 { 00228 c_idx[nSel]=SG_MALLOC(uint32_t, nz_dims); 00229 c_val[nSel]=SG_MALLOC(float64_t, nz_dims); 00230 00231 uint32_t idx=0; 00232 for(j=0; j < nDim; j++ ) 00233 { 00234 if(new_a[j] != 0) 00235 { 00236 c_idx[nSel][idx] = j; 00237 c_val[nSel][idx++] = new_a[j]; 00238 } 00239 } 00240 } 00241 00242 new_col_H[nSel] = sq_norm_a; 00243 00244 for(i=0; i < nSel; i++) 00245 { 00246 float64_t tmp = c_bias[nSel]*c_bias[i]; 00247 for(j=0; j < c_nzd[i]; j++) 00248 tmp += new_a[c_idx[i][j]]*c_val[i][j]; 00249 00250 new_col_H[i] = tmp; 00251 } 00252 //CMath::display_vector(new_col_H, nSel+1, "new_col_H"); 00253 //CMath::display_vector((int32_t*) c_idx[nSel], (int32_t) nz_dims, "c_idx"); 00254 //CMath::display_vector((float64_t*) c_val[nSel], nz_dims, "c_val"); 00255 return 0; 00256 } 00257 00258 int CSVMOcas::sort(float64_t* vals, float64_t* data, uint32_t size) 00259 { 00260 CMath::qsort_index(vals, data, size); 00261 return 0; 00262 } 00263 00264 /*---------------------------------------------------------------------- 00265 sparse_compute_output( output ) does the follwing: 00266 00267 output = data_X'*W; 00268 ----------------------------------------------------------------------*/ 00269 int CSVMOcas::compute_output(float64_t *output, void* ptr) 00270 { 00271 CSVMOcas* o = (CSVMOcas*) ptr; 00272 CDotFeatures* f=o->features; 00273 int32_t nData=f->get_num_vectors(); 00274 00275 float64_t* y = o->lab.vector; 00276 00277 f->dense_dot_range(output, 0, nData, y, o->w.vector, o->w.vlen, 0.0); 00278 00279 for (int32_t i=0; i<nData; i++) 00280 output[i]+=y[i]*o->bias; 00281 //CMath::display_vector(o->w, o->w.vlen, "w"); 00282 //CMath::display_vector(output, nData, "out"); 00283 return 0; 00284 } 00285 00286 /*---------------------------------------------------------------------- 00287 sq_norm_W = compute_W( alpha, nSel ) does the following: 00288 00289 oldW = W; 00290 W = sparse_A(:,1:nSel)'*alpha; 00291 sq_norm_W = W'*W; 00292 dp_WoldW = W'*oldW'; 00293 00294 ----------------------------------------------------------------------*/ 00295 void CSVMOcas::compute_W( 00296 float64_t *sq_norm_W, float64_t *dp_WoldW, float64_t *alpha, uint32_t nSel, 00297 void* ptr ) 00298 { 00299 CSVMOcas* o = (CSVMOcas*) ptr; 00300 uint32_t nDim= (uint32_t) o->w.vlen; 00301 CMath::swap(o->w.vector, o->old_w); 00302 float64_t* W=o->w.vector; 00303 float64_t* oldW=o->old_w; 00304 memset(W, 0, sizeof(float64_t)*nDim); 00305 float64_t old_bias=o->bias; 00306 float64_t bias=0; 00307 00308 float64_t** c_val = o->cp_value; 00309 uint32_t** c_idx = o->cp_index; 00310 uint32_t* c_nzd = o->cp_nz_dims; 00311 float64_t* c_bias = o->cp_bias; 00312 00313 for(uint32_t i=0; i<nSel; i++) 00314 { 00315 uint32_t nz_dims = c_nzd[i]; 00316 00317 if(nz_dims > 0 && alpha[i] > 0) 00318 { 00319 for(uint32_t j=0; j < nz_dims; j++) 00320 W[c_idx[i][j]] += alpha[i]*c_val[i][j]; 00321 } 00322 bias += c_bias[i]*alpha[i]; 00323 } 00324 00325 *sq_norm_W = SGVector<float64_t>::dot(W,W, nDim) + CMath::sq(bias); 00326 *dp_WoldW = SGVector<float64_t>::dot(W,oldW, nDim) + bias*old_bias; 00327 //SG_PRINT("nSel=%d sq_norm_W=%f dp_WoldW=%f\n", nSel, *sq_norm_W, *dp_WoldW); 00328 00329 o->bias = bias; 00330 o->old_bias = old_bias; 00331 } 00332 00333 void CSVMOcas::init() 00334 { 00335 use_bias=true; 00336 bufsize=3000; 00337 C1=1; 00338 C2=1; 00339 00340 epsilon=1e-3; 00341 method=SVM_OCAS; 00342 old_w=NULL; 00343 tmp_a_buf=NULL; 00344 cp_value=NULL; 00345 cp_index=NULL; 00346 cp_nz_dims=NULL; 00347 cp_bias=NULL; 00348 00349 primal_objective = 0.0; 00350 00351 m_parameters->add(&C1, "C1", "Cost constant 1."); 00352 m_parameters->add(&C2, "C2", "Cost constant 2."); 00353 m_parameters->add(&use_bias, "use_bias", 00354 "Indicates if bias is used."); 00355 m_parameters->add(&epsilon, "epsilon", "Convergence precision."); 00356 m_parameters->add(&bufsize, "bufsize", "Maximum number of cutting planes."); 00357 m_parameters->add((machine_int_t*) &method, "method", 00358 "SVMOcas solver type."); 00359 } 00360 00361 float64_t CSVMOcas::compute_primal_objective() const 00362 { 00363 return primal_objective; 00364 } 00365