00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 #include <shogun/lib/config.h>
00012 #ifdef HAVE_LAPACK
00013 #include <shogun/regression/svr/LibLinearRegression.h>
00014 #include <shogun/mathematics/Math.h>
00015 #include <shogun/labels/RegressionLabels.h>
00016 #include <shogun/optimization/liblinear/tron.h>
00017 #include <shogun/lib/Signal.h>
00018
00019 using namespace shogun;
00020
00021 CLibLinearRegression::CLibLinearRegression() :
00022 CLinearMachine()
00023 {
00024 init_defaults();
00025 }
00026
00027 CLibLinearRegression::CLibLinearRegression(float64_t C, CDotFeatures* feats, CLabels* labs) :
00028 CLinearMachine()
00029 {
00030 init_defaults();
00031 set_C(C);
00032 set_features(feats);
00033 set_labels(labs);
00034 }
00035
00036 void CLibLinearRegression::init_defaults()
00037 {
00038 set_C(1.0);
00039 set_epsilon(1e-2);
00040 set_max_iter(10000);
00041 set_use_bias(false);
00042 set_liblinear_regression_type(L2R_L1LOSS_SVR_DUAL);
00043 }
00044
00045 void CLibLinearRegression::register_parameters()
00046 {
00047 SG_ADD(&m_C, "m_C", "regularization constant",MS_AVAILABLE);
00048 SG_ADD(&m_epsilon, "m_epsilon", "tolerance epsilon",MS_NOT_AVAILABLE);
00049 SG_ADD(&m_epsilon, "m_tube_epsilon", "svr tube epsilon",MS_AVAILABLE);
00050 SG_ADD(&m_max_iter, "m_max_iter", "max number of iterations",MS_NOT_AVAILABLE);
00051 SG_ADD(&m_use_bias, "m_use_bias", "indicates whether bias should be used",MS_NOT_AVAILABLE);
00052 }
00053
00054 CLibLinearRegression::~CLibLinearRegression()
00055 {
00056 }
00057
00058 bool CLibLinearRegression::train_machine(CFeatures* data)
00059 {
00060 CSignal::clear_cancel();
00061
00062 if (data)
00063 set_features((CDotFeatures*)data);
00064
00065 ASSERT(features);
00066 ASSERT(m_labels && m_labels->get_label_type()==LT_REGRESSION);
00067
00068 int32_t num_train_labels=m_labels->get_num_labels();
00069 int32_t num_feat=features->get_dim_feature_space();
00070 int32_t num_vec=features->get_num_vectors();
00071
00072 if (num_vec!=num_train_labels)
00073 {
00074 SG_ERROR("number of vectors %d does not match "
00075 "number of training labels %d\n",
00076 num_vec, num_train_labels);
00077 }
00078
00079 if (m_use_bias)
00080 w=SGVector<float64_t>(SG_MALLOC(float64_t, num_feat+1), num_feat);
00081 else
00082 w=SGVector<float64_t>(num_feat);
00083
00084 problem prob;
00085 if (m_use_bias)
00086 {
00087 prob.n=w.vlen+1;
00088 memset(w.vector, 0, sizeof(float64_t)*(w.vlen+1));
00089 }
00090 else
00091 {
00092 prob.n=w.vlen;
00093 memset(w.vector, 0, sizeof(float64_t)*(w.vlen+0));
00094 }
00095 prob.l=num_vec;
00096 prob.x=features;
00097 prob.y=SG_MALLOC(float64_t, prob.l);
00098
00099 switch (m_liblinear_regression_type)
00100 {
00101 case L2R_L2LOSS_SVR:
00102 {
00103 double* Cs = SG_MALLOC(double, prob.l);
00104 for(int i = 0; i < prob.l; i++)
00105 Cs[i] = m_C;
00106
00107 function *fun_obj=new l2r_l2_svr_fun(&prob, Cs, m_tube_epsilon);
00108 CTron tron_obj(fun_obj, m_epsilon);
00109 tron_obj.tron(w.vector, m_max_train_time);
00110 delete fun_obj;
00111 SG_FREE(Cs);
00112 break;
00113
00114 }
00115 case L2R_L1LOSS_SVR_DUAL:
00116 solve_l2r_l1l2_svr(&prob);
00117 break;
00118 case L2R_L2LOSS_SVR_DUAL:
00119 solve_l2r_l1l2_svr(&prob);
00120 break;
00121 default:
00122 SG_ERROR("Error: unknown regression type\n");
00123 break;
00124 }
00125
00126 return true;
00127 }
00128
00129
00130
00131
00132
00133
00134
00135
00136
00137
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147
00148
00149
00150
00151
00152
00153 #undef GETI
00154 #define GETI(i) (0)
00155
00156
00157 void CLibLinearRegression::solve_l2r_l1l2_svr(const problem *prob)
00158 {
00159 int l = prob->l;
00160 double C = m_C;
00161 double p = m_tube_epsilon;
00162 int w_size = prob->n;
00163 double eps = m_epsilon;
00164 int i, s, iter = 0;
00165 int max_iter = 1000;
00166 int active_size = l;
00167 int *index = new int[l];
00168
00169 double d, G, H;
00170 double Gmax_old = CMath::INFTY;
00171 double Gmax_new, Gnorm1_new;
00172 double Gnorm1_init = 0.0;
00173 double *beta = new double[l];
00174 double *QD = new double[l];
00175 double *y = prob->y;
00176
00177
00178 double lambda[1], upper_bound[1];
00179 lambda[0] = 0.5/C;
00180 upper_bound[0] = CMath::INFTY;
00181
00182 if(m_liblinear_regression_type == L2R_L1LOSS_SVR_DUAL)
00183 {
00184 lambda[0] = 0;
00185 upper_bound[0] = C;
00186 }
00187
00188
00189
00190 for(i=0; i<l; i++)
00191 beta[i] = 0;
00192
00193 for(i=0; i<w_size; i++)
00194 w[i] = 0;
00195 for(i=0; i<l; i++)
00196 {
00197 QD[i] = prob->x->dot(i, prob->x,i);
00198 prob->x->add_to_dense_vec(beta[i], i, w.vector, w_size);
00199
00200 if (prob->use_bias)
00201 w.vector[w_size]+=beta[i];
00202
00203 index[i] = i;
00204 }
00205
00206
00207 while(iter < max_iter)
00208 {
00209 Gmax_new = 0;
00210 Gnorm1_new = 0;
00211
00212 for(i=0; i<active_size; i++)
00213 {
00214 int j = i+rand()%(active_size-i);
00215 CMath::swap(index[i], index[j]);
00216 }
00217
00218 for(s=0; s<active_size; s++)
00219 {
00220 i = index[s];
00221 G = -y[i] + lambda[GETI(i)]*beta[i];
00222 H = QD[i] + lambda[GETI(i)];
00223
00224 G += prob->x->dense_dot(i, w.vector, w_size);
00225 if (prob->use_bias)
00226 G+=w.vector[w_size];
00227
00228 double Gp = G+p;
00229 double Gn = G-p;
00230 double violation = 0;
00231 if(beta[i] == 0)
00232 {
00233 if(Gp < 0)
00234 violation = -Gp;
00235 else if(Gn > 0)
00236 violation = Gn;
00237 else if(Gp>Gmax_old && Gn<-Gmax_old)
00238 {
00239 active_size--;
00240 CMath::swap(index[s], index[active_size]);
00241 s--;
00242 continue;
00243 }
00244 }
00245 else if(beta[i] >= upper_bound[GETI(i)])
00246 {
00247 if(Gp > 0)
00248 violation = Gp;
00249 else if(Gp < -Gmax_old)
00250 {
00251 active_size--;
00252 CMath::swap(index[s], index[active_size]);
00253 s--;
00254 continue;
00255 }
00256 }
00257 else if(beta[i] <= -upper_bound[GETI(i)])
00258 {
00259 if(Gn < 0)
00260 violation = -Gn;
00261 else if(Gn > Gmax_old)
00262 {
00263 active_size--;
00264 CMath::swap(index[s], index[active_size]);
00265 s--;
00266 continue;
00267 }
00268 }
00269 else if(beta[i] > 0)
00270 violation = fabs(Gp);
00271 else
00272 violation = fabs(Gn);
00273
00274 Gmax_new = CMath::max(Gmax_new, violation);
00275 Gnorm1_new += violation;
00276
00277
00278 if(Gp < H*beta[i])
00279 d = -Gp/H;
00280 else if(Gn > H*beta[i])
00281 d = -Gn/H;
00282 else
00283 d = -beta[i];
00284
00285 if(fabs(d) < 1.0e-12)
00286 continue;
00287
00288 double beta_old = beta[i];
00289 beta[i] = CMath::min(CMath::max(beta[i]+d, -upper_bound[GETI(i)]), upper_bound[GETI(i)]);
00290 d = beta[i]-beta_old;
00291
00292 if(d != 0)
00293 {
00294 prob->x->add_to_dense_vec(d, i, w.vector, w_size);
00295
00296 if (prob->use_bias)
00297 w.vector[w_size]+=d;
00298 }
00299 }
00300
00301 if(iter == 0)
00302 Gnorm1_init = Gnorm1_new;
00303 iter++;
00304
00305 SG_SABS_PROGRESS(Gnorm1_new, -CMath::log10(Gnorm1_new), -CMath::log10(eps*Gnorm1_init), -CMath::log10(Gnorm1_init), 6);
00306
00307 if(Gnorm1_new <= eps*Gnorm1_init)
00308 {
00309 if(active_size == l)
00310 break;
00311 else
00312 {
00313 active_size = l;
00314 Gmax_old = CMath::INFTY;
00315 continue;
00316 }
00317 }
00318
00319 Gmax_old = Gmax_new;
00320 }
00321
00322 SG_DONE();
00323 SG_INFO("\noptimization finished, #iter = %d\n", iter);
00324 if(iter >= max_iter)
00325 SG_INFO("\nWARNING: reaching max number of iterations\nUsing -s 11 may be faster\n\n");
00326
00327
00328 double v = 0;
00329 int nSV = 0;
00330 for(i=0; i<w_size; i++)
00331 v += w[i]*w[i];
00332 v = 0.5*v;
00333 for(i=0; i<l; i++)
00334 {
00335 v += p*fabs(beta[i]) - y[i]*beta[i] + 0.5*lambda[GETI(i)]*beta[i]*beta[i];
00336 if(beta[i] != 0)
00337 nSV++;
00338 }
00339
00340 SG_INFO("Objective value = %lf\n", v);
00341 SG_INFO("nSV = %d\n",nSV);
00342
00343 delete [] beta;
00344 delete [] QD;
00345 delete [] index;
00346 }
00347
00348 #endif