Mercurial > hg > camir-aes2014

annotate toolboxes/SVM-light/src/svm_loqo.c @ 0:e9a9cd732c1e tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
first hg version after svn
author wolffd
date Tue, 10 Feb 2015 15:05:51 +0000
parents
children
rev   line source
wolffd@0 1 /***********************************************************************/
wolffd@0 2 /* */
wolffd@0 3 /* svm_loqo.c */
wolffd@0 4 /* */
wolffd@0 5 /* Interface to the PR_LOQO optimization package for SVM. */
wolffd@0 6 /* */
wolffd@0 7 /* Author: Thorsten Joachims */
wolffd@0 8 /* Date: 19.07.99 */
wolffd@0 9 /* */
wolffd@0 10 /* Copyright (c) 1999 Universitaet Dortmund - All rights reserved */
wolffd@0 11 /* */
wolffd@0 12 /* This software is available for non-commercial use only. It must */
wolffd@0 13 /* not be modified and distributed without prior permission of the */
wolffd@0 14 /* author. The author is not responsible for implications from the */
wolffd@0 15 /* use of this software. */
wolffd@0 16 /* */
wolffd@0 17 /***********************************************************************/
wolffd@0 18
wolffd@0 19 # include <math.h>
wolffd@0 20 # include "pr_loqo/pr_loqo.h"
wolffd@0 21 # include "svm_common.h"
wolffd@0 22
wolffd@0 23 /* Common Block Declarations */
wolffd@0 24
wolffd@0 25 long verbosity;
wolffd@0 26
wolffd@0 27 /* /////////////////////////////////////////////////////////////// */
wolffd@0 28
wolffd@0 29 # define DEF_PRECISION_LINEAR 1E-8
wolffd@0 30 # define DEF_PRECISION_NONLINEAR 1E-14
wolffd@0 31
wolffd@0 32 double *optimize_qp();
wolffd@0 33 double *primal=0,*dual=0;
wolffd@0 34 double init_margin=0.15;
wolffd@0 35 long init_iter=500,precision_violations=0;
wolffd@0 36 double model_b;
wolffd@0 37 double opt_precision=DEF_PRECISION_LINEAR;
wolffd@0 38
wolffd@0 39 /* /////////////////////////////////////////////////////////////// */
wolffd@0 40
wolffd@0 41 void *my_malloc();
wolffd@0 42
wolffd@0 43 double *optimize_qp(qp,epsilon_crit,nx,threshold,learn_parm)
wolffd@0 44 QP *qp;
wolffd@0 45 double *epsilon_crit;
wolffd@0 46 long nx; /* Maximum number of variables in QP */
wolffd@0 47 double *threshold;
wolffd@0 48 LEARN_PARM *learn_parm;
wolffd@0 49 /* start the optimizer and return the optimal values */
wolffd@0 50 {
wolffd@0 51 register long i,j,result;
wolffd@0 52 double margin,obj_before,obj_after;
wolffd@0 53 double sigdig,dist,epsilon_loqo;
wolffd@0 54 int iter;
wolffd@0 55
wolffd@0 56 if(!primal) { /* allocate memory at first call */
wolffd@0 57 primal=(double *)my_malloc(sizeof(double)*nx*3);
wolffd@0 58 dual=(double *)my_malloc(sizeof(double)*(nx*2+1));
wolffd@0 59 }
wolffd@0 60
wolffd@0 61 if(verbosity>=4) { /* really verbose */
wolffd@0 62 printf("\n\n");
wolffd@0 63 for(i=0;i<qp->opt_n;i++) {
wolffd@0 64 printf("%f: ",qp->opt_g0[i]);
wolffd@0 65 for(j=0;j<qp->opt_n;j++) {
wolffd@0 66 printf("%f ",qp->opt_g[i*qp->opt_n+j]);
wolffd@0 67 }
wolffd@0 68 printf(": a%ld=%.10f < %f",i,qp->opt_xinit[i],qp->opt_up[i]);
wolffd@0 69 printf(": y=%f\n",qp->opt_ce[i]);
wolffd@0 70 }
wolffd@0 71 for(j=0;j<qp->opt_m;j++) {
wolffd@0 72 printf("EQ-%ld: %f*a0",j,qp->opt_ce[j]);
wolffd@0 73 for(i=1;i<qp->opt_n;i++) {
wolffd@0 74 printf(" + %f*a%ld",qp->opt_ce[i],i);
wolffd@0 75 }
wolffd@0 76 printf(" = %f\n\n",-qp->opt_ce0[0]);
wolffd@0 77 }
wolffd@0 78 }
wolffd@0 79
wolffd@0 80 obj_before=0; /* calculate objective before optimization */
wolffd@0 81 for(i=0;i<qp->opt_n;i++) {
wolffd@0 82 obj_before+=(qp->opt_g0[i]*qp->opt_xinit[i]);
wolffd@0 83 obj_before+=(0.5*qp->opt_xinit[i]*qp->opt_xinit[i]*qp->opt_g[i*qp->opt_n+i]);
wolffd@0 84 for(j=0;j<i;j++) {
wolffd@0 85 obj_before+=(qp->opt_xinit[j]*qp->opt_xinit[i]*qp->opt_g[j*qp->opt_n+i]);
wolffd@0 86 }
wolffd@0 87 }
wolffd@0 88
wolffd@0 89 result=STILL_RUNNING;
wolffd@0 90 qp->opt_ce0[0]*=(-1.0);
wolffd@0 91 /* Run pr_loqo. If a run fails, try again with parameters which lead */
wolffd@0 92 /* to a slower, but more robust setting. */
wolffd@0 93 for(margin=init_margin,iter=init_iter;
wolffd@0 94 (margin<=0.9999999) && (result!=OPTIMAL_SOLUTION);) {
wolffd@0 95 sigdig=-log10(opt_precision);
wolffd@0 96
wolffd@0 97 result=pr_loqo((int)qp->opt_n,(int)qp->opt_m,
wolffd@0 98 (double *)qp->opt_g0,(double *)qp->opt_g,
wolffd@0 99 (double *)qp->opt_ce,(double *)qp->opt_ce0,
wolffd@0 100 (double *)qp->opt_low,(double *)qp->opt_up,
wolffd@0 101 (double *)primal,(double *)dual,
wolffd@0 102 (int)(verbosity-2),
wolffd@0 103 (double)sigdig,(int)iter,
wolffd@0 104 (double)margin,(double)(qp->opt_up[0])/4.0,(int)0);
wolffd@0 105
wolffd@0 106 if(isnan(dual[0])) { /* check for choldc problem */
wolffd@0 107 if(verbosity>=2) {
wolffd@0 108 printf("NOTICE: Restarting PR_LOQO with more conservative parameters.\n");
wolffd@0 109 }
wolffd@0 110 if(init_margin<0.80) { /* become more conservative in general */
wolffd@0 111 init_margin=(4.0*margin+1.0)/5.0;
wolffd@0 112 }
wolffd@0 113 margin=(margin+1.0)/2.0;
wolffd@0 114 (opt_precision)*=10.0; /* reduce precision */
wolffd@0 115 if(verbosity>=2) {
wolffd@0 116 printf("NOTICE: Reducing precision of PR_LOQO.\n");
wolffd@0 117 }
wolffd@0 118 }
wolffd@0 119 else if(result!=OPTIMAL_SOLUTION) {
wolffd@0 120 iter+=2000;
wolffd@0 121 init_iter+=10;
wolffd@0 122 (opt_precision)*=10.0; /* reduce precision */
wolffd@0 123 if(verbosity>=2) {
wolffd@0 124 printf("NOTICE: Reducing precision of PR_LOQO due to (%ld).\n",result);
wolffd@0 125 }
wolffd@0 126 }
wolffd@0 127 }
wolffd@0 128
wolffd@0 129 if(qp->opt_m) /* Thanks to Alex Smola for this hint */
wolffd@0 130 model_b=dual[0];
wolffd@0 131 else
wolffd@0 132 model_b=0;
wolffd@0 133
wolffd@0 134 /* Check the precision of the alphas. If results of current optimization */
wolffd@0 135 /* violate KT-Conditions, relax the epsilon on the bounds on alphas. */
wolffd@0 136 epsilon_loqo=1E-10;
wolffd@0 137 for(i=0;i<qp->opt_n;i++) {
wolffd@0 138 dist=-model_b*qp->opt_ce[i];
wolffd@0 139 dist+=(qp->opt_g0[i]+1.0);
wolffd@0 140 for(j=0;j<i;j++) {
wolffd@0 141 dist+=(primal[j]*qp->opt_g[j*qp->opt_n+i]);
wolffd@0 142 }
wolffd@0 143 for(j=i;j<qp->opt_n;j++) {
wolffd@0 144 dist+=(primal[j]*qp->opt_g[i*qp->opt_n+j]);
wolffd@0 145 }
wolffd@0 146 /* printf("LOQO: a[%d]=%f, dist=%f, b=%f\n",i,primal[i],dist,dual[0]); */
wolffd@0 147 if((primal[i]<(qp->opt_up[i]-epsilon_loqo)) && (dist < (1.0-(*epsilon_crit)))) {
wolffd@0 148 epsilon_loqo=(qp->opt_up[i]-primal[i])*2.0;
wolffd@0 149 }
wolffd@0 150 else if((primal[i]>(0+epsilon_loqo)) && (dist > (1.0+(*epsilon_crit)))) {
wolffd@0 151 epsilon_loqo=primal[i]*2.0;
wolffd@0 152 }
wolffd@0 153 }
wolffd@0 154
wolffd@0 155 for(i=0;i<qp->opt_n;i++) { /* clip alphas to bounds */
wolffd@0 156 if(primal[i]<=(0+epsilon_loqo)) {
wolffd@0 157 primal[i]=0;
wolffd@0 158 }
wolffd@0 159 else if(primal[i]>=(qp->opt_up[i]-epsilon_loqo)) {
wolffd@0 160 primal[i]=qp->opt_up[i];
wolffd@0 161 }
wolffd@0 162 }
wolffd@0 163
wolffd@0 164 obj_after=0; /* calculate objective after optimization */
wolffd@0 165 for(i=0;i<qp->opt_n;i++) {
wolffd@0 166 obj_after+=(qp->opt_g0[i]*primal[i]);
wolffd@0 167 obj_after+=(0.5*primal[i]*primal[i]*qp->opt_g[i*qp->opt_n+i]);
wolffd@0 168 for(j=0;j<i;j++) {
wolffd@0 169 obj_after+=(primal[j]*primal[i]*qp->opt_g[j*qp->opt_n+i]);
wolffd@0 170 }
wolffd@0 171 }
wolffd@0 172
wolffd@0 173 /* if optimizer returned NAN values, reset and retry with smaller */
wolffd@0 174 /* working set. */
wolffd@0 175 if(isnan(obj_after) || isnan(model_b)) {
wolffd@0 176 for(i=0;i<qp->opt_n;i++) {
wolffd@0 177 primal[i]=qp->opt_xinit[i];
wolffd@0 178 }
wolffd@0 179 model_b=0;
wolffd@0 180 if(learn_parm->svm_maxqpsize>2) {
wolffd@0 181 learn_parm->svm_maxqpsize--; /* decrease size of qp-subproblems */
wolffd@0 182 }
wolffd@0 183 }
wolffd@0 184
wolffd@0 185 if(obj_after >= obj_before) { /* check whether there was progress */
wolffd@0 186 (opt_precision)/=100.0;
wolffd@0 187 precision_violations++;
wolffd@0 188 if(verbosity>=2) {
wolffd@0 189 printf("NOTICE: Increasing Precision of PR_LOQO.\n");
wolffd@0 190 }
wolffd@0 191 }
wolffd@0 192
wolffd@0 193 if(precision_violations > 500) {
wolffd@0 194 (*epsilon_crit)*=10.0;
wolffd@0 195 precision_violations=0;
wolffd@0 196 if(verbosity>=1) {
wolffd@0 197 printf("\nWARNING: Relaxing epsilon on KT-Conditions.\n");
wolffd@0 198 }
wolffd@0 199 }
wolffd@0 200
wolffd@0 201 (*threshold)=model_b;
wolffd@0 202
wolffd@0 203 if(result!=OPTIMAL_SOLUTION) {
wolffd@0 204 printf("\nERROR: PR_LOQO did not converge. \n");
wolffd@0 205 return(qp->opt_xinit);
wolffd@0 206 }
wolffd@0 207 else {
wolffd@0 208 return(primal);
wolffd@0 209 }
wolffd@0 210 }
wolffd@0 211