/*
File name: p_svm.c
Created by: Ljubomir Buturovic
Created: 08/23/2002
Purpose: SVM model selection, learning and prediction, using the
LIBSVM library.
*/
/*
Copyright 2004 Ljubomir J. Buturovic
Permission is hereby granted, free of charge, to any person
obtaining a copy of this software and associated documentation files
(the "Software"), to deal in the Software without restriction,
including without limitation the rights to use, copy, modify, merge,
publish, distribute, sublicense, and/or sell copies of the Software,
and to permit persons to whom the Software is furnished to do so,
subject to the following conditions:
The above copyright notice and this permission notice shall be
included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS
BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN
ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
SOFTWARE.
*/
static char rcsid[] = "$Id: p_svm.c,v 1.127 2006/05/24 06:06:06 ljubomir Exp $";
#include <stdlib.h>
#include <errno.h>
#include <unistd.h>
#include <string.h>
#include <float.h>
#include <math.h>
#include "simplex.h"
#include "lau.h"
#include "lmat.h"
#include "pcl_svm.h"
#include "pau.h"
#include "svm.h"
#include "bagging.h"
#include "xlearn.h"
#include "xpar.h"
/*
This function collects SVM learning parameters from the user, and
then calls svm_train() which performs the SVM training.
In case of successful completion, set 'errc' to 0. In case of error,
set errc to 'errno'. If error is file access error, set 'file'.
*/
void p_svm_learn(int *errc, char **xname, int dbg)
{
int status;
char *fname;
FILE *outdev;
FILE *fdbg = (FILE *) 0;
struct svm_problem *problem;
struct svm_parameter *parameters;
struct svm_model *model;
*errc = 0;
fname = (char *) 0; /* avoid compiler warning */
outdev = stdout;
fflush(outdev);
if (dbg)
fdbg = fopen(PCP_DBG, "a");
problem = create_problem(tds);
if (problem)
{
parameters = calloc(1, sizeof(struct svm_parameter));
if (parameters)
{
status = input_svm(outdev, tds->c, tds->nv, problem, parameters,
(char **) 0, (int *) 0, 0, 1, PCP_SVM_K_NONE,
(int *) 0);
parameters->probability = 1;
if (!status)
{
fname = input_filename(PCP_UMSG_SVM_FNAME, PCP_SVM, outdev);
inverse_video();
model = svm_train(problem, parameters);
if (model == (struct svm_model *) 0)
status = -1;
else
status = svm_save_model(fname, model);
if (status == -1)
{
*errc = errno;
*xname = strdup(fname);
}
}
else
{
*errc = errno;
*xname = strdup(fname);
}
free(parameters);
}
else
*errc = errno;
}
else
*errc = errno;
reset_video();
}
#define PCP_NU_SVC 1
#define PCP_C_SVC 2
/*
Kernel types are internally coded starting from 1 (because 0 is
reserved for the default kernel type). This function converts the
internal kernel types to LIBSVM kernel types.
*/
int convert_kernel_type(int type)
{
int kernel_type;
kernel_type = PCP_SVM_K_RBF;
if (type == PCP_SVM_K_LINEAR)
kernel_type = LINEAR;
else if (type == PCP_SVM_K_POLY)
kernel_type = POLY;
else if (type == PCP_SVM_K_RBF)
kernel_type = RBF;
else if (type == PCP_SVM_K_SIGMOID)
kernel_type = SIGMOID;
return kernel_type;
}
int input_svm_type(char *msg, FILE *outdev)
{
int min_range;
int max_range;
int type;
min_range = PCP_NU_SVC;
max_range = PCP_C_SVC;
sprintf(msg, PCP_UMSG_SVM_TYPE, PCP_NU_SVC, PCP_C_SVC, PCP_C_SVC);
type = input_integer(stdin, outdev, msg, PCP_QLEN, &max_range, &min_range, &max_range);
if (type == PCP_NU_SVC)
type = NU_SVC;
else if (type == PCP_C_SVC)
type = C_SVC;
return type;
}
/*
Accept kernel type from stdin.
*/
int input_kernel_type(char *msg, FILE *outdev)
{
int kernel_type;
int min_range;
int max_range;
int dflt;
min_range = PCP_SVM_K_LINEAR;
max_range = PCP_SVM_K_SIGMOID;
dflt = PCP_SVM_K_RBF;
sprintf(msg, PCP_UMSG_KERNEL_TYPE, PCP_SVM_K_LINEAR, PCP_SVM_K_POLY, PCP_SVM_K_RBF,
PCP_SVM_K_SIGMOID, dflt);
kernel_type = input_integer(stdin, outdev, msg, PCP_QLEN, &dflt, &min_range, &max_range);
return kernel_type;
}
/*
Accept kernel type from a limited selection (only the types
supported by the model selection functions pcp_svm_grid() and
pcp_svm_simplex()).
*/
int input_kernel(char *msg, FILE *outdev)
{
int kernel_type;
int dflt;
int choices[2];
int cntr;
cntr = 0;
choices[cntr++] = PCP_SVM_K_LINEAR;
choices[cntr++] = PCP_SVM_K_RBF;
dflt = PCP_SVM_K_RBF;
sprintf(msg, PCP_UMSG_KERNEL, PCP_SVM_K_LINEAR, PCP_SVM_K_RBF, dflt);
kernel_type = input_choice(stdin, outdev, msg, &dflt, choices, cntr);
return kernel_type;
}
void input_class_costs(struct svm_parameter *parameters, int c)
{
int min_range;
int max_range;
int dflt;
int ixd;
int i;
float fmin;
float fdflt;
float fxd;
char *msg;
min_range = 0;
max_range = 1;
dflt = 0;
ixd = input_integer(stdin, stdout, PCP_UMSG_CCOSTS, PCP_QLEN, &dflt, &min_range, &max_range);
if (ixd == 1)
{
parameters->nr_weight = c;
parameters->weight_label = malloc(c*sizeof(int));
parameters->weight = malloc(c*sizeof(double));
msg = malloc((PCP_QLEN+1)*sizeof(char));
for (i = 0; i < c; i++)
{
fmin = 0.0;
fdflt = 1.0;
sprintf(msg, PCP_UMSG_CLASS_COST, i+1, fdflt);
fxd = input_float(stdin, stdout, msg, PCP_QLEN, &fdflt, &fmin, (float *) 0);
parameters->weight_label[i] = i+1;
parameters->weight[i] = fxd;
}
free(msg);
}
}
/*
Set the SVM parameters we decided to keep as defaults.
*/
void svm_param_defaults(void *svm_parameters)
{
float fdflt;
struct svm_parameter *svm_params;
svm_params = (struct svm_parameter *) svm_parameters;
svm_params->probability = PCP_SVM_DFLT_PROB;
svm_params->cache_size = PCP_SVM_DFLT_CACHE;
if (svm_params->svm_type == NU_SVC) /* default values per LIBSVM README file */
fdflt = PCP_SVM_DFLT_EPS_NU;
else
fdflt = PCP_SVM_DFLT_EPS_C;
svm_params->eps = fdflt; /* decided to always use default values. ljb, 03/06/2005 */
svm_params->shrinking = 1; /* Decided to always use the shrinking heuristics */
}
/*
Compute max. feasible value of nu for NU-SVM. Use Proposition 3 in
Section 6 of
http://www.csie.ntu.edu.tw/~cjlin/papers/nusvmtutorial.pdf.
*/
float get_svm_nu_max(int *nd, int c)
{
int i;
int j;
int n1;
int n2;
int nx;
float nu;
float numax;
numax = 1.0;
if (nd)
{
for (i = 0; i < c-1; i++)
for (j = i+1; j < c; j++)
{
n1 = nd[i];
n2 = nd[j];
nx = n1;
if (n2 < n1)
nx = n2;
nu = 2.0*nx/(n1+n2);
if (nu < numax)
numax = nu;
}
}
return numax;
}
/*
Compute max. feasible value of nu for `nxval' NU-SVM
cross-validation experiments.
*/
float get_svm_nu_xmax(int *nd, int c, int nxval)
{
int i;
int j;
int n1;
int n2;
float nu;
float numax;
int ndx[2];
numax = 1.0;
if (nd)
{
for (i = 0; i < c-1; i++)
for (j = i+1; j < c; j++)
{
n1 = nd[i];
n2 = nd[j];
/*
These steps go through all four possibilities of
learning subset sizes for classes i and j, and pick the
lowest numax of the four possibilities.
*/
ndx[0] = n1-n1/nxval;
ndx[1] = n2-n2/nxval;
nu = get_svm_nu_max(ndx, 2);
if (nu < numax)
numax = nu;
ndx[0] -= 1;
nu = get_svm_nu_max(ndx, 2);
if (nu < numax)
numax = nu;
ndx[0] += 1;
ndx[1] -= 1;
nu = get_svm_nu_max(ndx, 2);
if (nu < numax)
numax = nu;
ndx[0] -= 1;
nu = get_svm_nu_max(ndx, 2);
if (nu < numax)
numax = nu;
}
}
return numax;
}
/*
Read SVM parameters from stdin. Check consistency with 'problem'.
If `full_par' is not 0, read C/nu, gamma, degree and coef0
parameters, depending on the SVM type. Otherwise, if `ktype' is
PCP_SVM_K_NONE, read only the kernel type. Otherwise, we don't need
these parameter values (presumably because we are called from a
function which will optimize them).
If `method' is not NULL, read the ensemble variant of PALG_SVM
(i.e., PALG_BAG_SVM or PALG_ADABOOST_SVM).
*/
int input_svm(FILE *outdev, int c, int nvec, struct svm_problem *problem,
struct svm_parameter *parameters, char **fname, int *nxval,
int max_nxval, int full_par, int ktype, int *ensemble_method)
{
int i;
int status = 0;
int ixd;
int kernel_type;
int min_range;
int max_range;
int done;
float fxd;
float fmin;
float fmax;
float fdflt;
char *pcheck;
char *msg;
FILE *fptr;
clear_screen();
cursor_on();
if (ensemble_method)
*ensemble_method = input_ensemble_method(PALG_SVM);
msg = malloc((PCP_QLEN+1)*sizeof(char));
done = 0;
while (done == 0)
{
parameters->svm_type = input_svm_type(msg, outdev);
if (full_par)
{
kernel_type = input_kernel_type(msg, outdev);
ixd = convert_kernel_type(kernel_type);
parameters->kernel_type = ixd;
}
else if (ktype != PCP_SVM_K_NONE)
{
kernel_type = input_kernel(msg, outdev);
parameters->kernel_type = convert_kernel_type(kernel_type);
}
else
{
kernel_type = PCP_SVM_K_RBF;
parameters->kernel_type = convert_kernel_type(kernel_type);
}
if (c > 2)
{
min_range = 0;
max_range = 1;
sprintf(msg, PCP_UMSG_INDET, max_range, min_range, min_range);
/* TBD-indet: disable indeterminate-region multi-class problem algorithm for now */
/*parameters->indet = input_integer(stdin, outdev, msg, PCP_QLEN, &dflt, &min_range, &max_range);*/
}
if (parameters->svm_type == NU_SVC && full_par)
{
fmin = 0.0;
fmax = get_svm_nu_max(tds->nd, tds->c);
if (PCP_SVM_DFLT_NU < fmax)
fdflt = PCP_SVM_DFLT_NU;
else
fdflt = 0.5*fdflt;
sprintf(msg, "%s (%4.2f..%4.2f) [%4.2f]:", PCP_UMSG_NU, fmin, fmax, fdflt);
fxd = input_float(stdin, outdev, msg, PCP_QLEN, &fdflt, &fmin, &fmax);
parameters->nu = fxd;
}
if (parameters->svm_type == C_SVC && full_par)
{
fdflt = PCP_SVM_DFLT_C;
fmin = 0.0;
sprintf(msg, PCP_UMSG_COST, fdflt);
fxd = input_float(stdin, outdev, msg, PCP_QLEN, &fdflt, &fmin, (float *) 0);
parameters->C = fxd;
}
if (parameters->svm_type == C_SVC)
input_class_costs(parameters, c);
if (full_par)
{
if (parameters->kernel_type == POLY)
{
fdflt = PCP_SVM_DFLT_DEGREE;
fmin = 1.0;
sprintf(msg, PCP_UMSG_DEGREE, fdflt);
fxd = input_float(stdin, outdev, msg, PCP_QLEN, &fdflt, &fmin, (float *) 0);
parameters->degree = fxd;
}
if ((parameters->kernel_type == POLY) || (parameters->kernel_type == RBF) ||
(parameters->kernel_type == SIGMOID))
{
/*fdflt = 1.0/(tds->d);*/
fdflt = PCP_SVM_DFLT_GAMMA;
sprintf(msg, PCP_UMSG_GAMMA, fdflt);
fxd = input_float(stdin, outdev, msg, PCP_QLEN, &fdflt, (float *) 0, (float *) 0);
parameters->gamma = fxd;
}
if ((parameters->kernel_type == POLY) || (parameters->kernel_type == SIGMOID))
{
fdflt = PCP_SVM_DFLT_COEF0;
sprintf(msg, PCP_UMSG_COEF0, fdflt);
fxd = input_float(stdin, outdev, msg, PCP_QLEN, &fdflt, (float *) 0, (float *) 0);
parameters->coef0 = fxd;
}
}
svm_param_defaults(parameters);
pcheck = (char *) 0;
if (problem)
pcheck = (char *) svm_check_parameter(problem, parameters);
if (pcheck)
{
sprintf(msg, PCP_UMSG_BAD_SVM_PARAMS, pcheck);
print_line(outdev, msg, PCP_QLEN);
fprintf(outdev, "\n");
}
else
done = 1;
}
if (nxval)
*nxval = input_nxval(stdin, outdev, max_nxval);
if (fname)
{
if (ensemble_method)
*fname = input_filename(PCP_UMSG_SVM_FNAME, PCP_SVM, outdev);
else
*fname = input_filename(PCP_UMSG_XSV, PCP_XSV, outdev);
}
if (fname)
{
fptr = fopen(*fname, "w");
if (fptr)
{
fclose(fptr);
unlink(*fname);
}
else
{
i = errno; /* for debugging */
status = -1;
}
}
vx_free(msg);
return status;
}
/*
Save training data set in LIBSVM data format. See LIBSVM
documentation for the format description.
In case of success, set 'errc' to 0, otherwise set it to 'errno'. If
error is file access error, store the name of the offending file in
'xname'.
*/
void p_svm_save(int *errc, char **xname)
{
int i;
int j;
int type;
int min_range;
int max_range;
char *fname;
char *msg;
struct svm_problem *problem;
struct svm_node *node;
FILE *fptr;
*errc = 0;
clear_screen();
cursor_on();
msg = malloc((PCP_QLEN+1)*sizeof(char));
min_range = 0;
max_range = 1;
sprintf(msg, " Convert training (%d) or test (%d) dataset [%d]:",
min_range, max_range, min_range);
type = input_integer(stdin, stdout, msg, PCP_QLEN, &min_range, &min_range,
&max_range);
free(msg);
fname = input_filename(PCP_UMSG_ONAME, PCP_LVM, stdout);
if (type == 0)
problem = create_problem(tds);
else
problem = create_problem(teds);
if (problem)
{
fptr = fopen(fname, "w");
if (fptr)
{
for (i = 0; i < problem->l; i++)
{
fprintf(fptr, "%d ", (int) problem->y[i]);
node = problem->x[i];
j = 0;
while (node[j].index != -1)
{
fprintf(fptr, "%d:%g ", j+1, node[j].value);
j++;
}
fprintf(fptr, "\n");
}
fclose(fptr);
}
else
{
*errc = errno;
*xname = strdup(fname);
}
}
else
*errc = errno;
}
/*
Classify test data set using SVM classifier stored in 'svm_fname'
(the file name provided by user), and display classification
results.
In case of successful completion, set 'errc' to 0. Otherwise, set
errc to 'errno'. In case of file error, store the name of the file
causing the error in 'xname'.
*/
void p_svm_predict(int *errc, char **xname)
{
int i;
int j;
int verbose;
int min_range;
int max_range;
int predicted_class;
int nmodels;
int prob_flag;
float *predictions;
char *svm_fname;
char *output_fname;
float *weights;
struct svm_model *model;
struct svm_node *svm_vector;
struct svm_model **models;
FILE *outdev;
FILE *svm_fptr;
FILE *output_fptr;
*errc = 0;
verbose = -1;
clear_screen();
outdev = stdout;
cursor_on();
svm_fname = input_filename(PCP_UMSG_SVM_FNAME, PCP_SVM, outdev);
svm_fptr = fopen(svm_fname, "r");
if (svm_fptr)
{
fclose(svm_fptr);
output_fname = strdup(PCP_RCL);
output_fptr = fopen(output_fname, "w");
if (output_fptr)
{
models = pcl_svm_load(svm_fname, &nmodels, &weights);
if (models)
{
predictions = calloc(tds->c, sizeof(float));
teds->prediction = malloc(teds->nv*sizeof(int));
teds->ps = dmx_alloc(teds->nv, teds->c);
for (i = 0; i < teds->nv; i++)
{
svm_vector = create_svm_vector(teds->x[i], teds->d);
fvec_set(predictions, tds->c, 0.0);
for (j = 0; j < nmodels; j++)
{
model = models[j];
/*
Check if the model contains probability
information.
*/
prob_flag = svm_check_probability_model(model);
if (prob_flag)
predicted_class = svm_predict_probability(model, svm_vector,
teds->ps[i])-1;
else
predicted_class = svm_predict(model, svm_vector)-1;
predictions[predicted_class] += weights[j];
}
predicted_class = fvec_argmax(predictions, tds->c);
free(svm_vector);
teds->prediction[i] = predicted_class;
write_rcl(output_fptr, i, teds, tds);
}
free(predictions);
for (j = 0; j < nmodels; j++)
svm_destroy_model(models[j]);
free(models);
fclose(output_fptr);
if (teds->c == tds->c)
{
min_range = 0;
max_range = 1;
verbose = input_integer(stdin, outdev, OUTPUT_MSG, PCP_QLEN,
&min_range, &min_range, &max_range);
}
/*
Display results.
*/
predict_disp(teds, verbose, PALG_SVM);
pwait();
}
else
{
*errc = PERR_UNRECOGNIZED_SVM;
*xname = strdup(svm_fname);
fclose(output_fptr);
unlink(output_fname);
}
}
else
{
*errc = errno;
*xname = strdup(output_fname);
}
}
else
{
*errc = errno;
*xname = strdup(svm_fname);
}
}
#define PCP_GAMMA_LEFT -8
#define PCP_GAMMA_RIGHT 0
#define PCP_GAMMA_MIDDLE -2
#define PCP_C_LEFT -1
#define PCP_C_RIGHT 4
#define PCP_C_MIDDLE 1.5
#define PCP_NU_LEFT 0.01
#define PCP_NU_FACTOR 0.99 /* numeric safety factor for NU-SVM feasible region */
#define PCP_NSTEP 10 /* number of steps for grid search */
#define PCP_RSTEP 5 /* number of steps for refinement of the grid search */
#define PCP_RFACTOR 0.05 /* refinement factor */
/*#define PCP_GAMMA_LEFT -6.75*/ /* for testing */
/*#define PCP_C_LEFT 3.609375*/ /* for testing */
/*
Transplant optimal parameters to xpar_crit_parameters->options.
*/
static void set_optimal(struct xpar_crit_parameters *xpar_parameters)
{
if (((struct svm_parameter *) xpar_parameters->options)->svm_type == C_SVC)
((struct svm_parameter *) xpar_parameters->options)->C = pow(10, xpar_parameters->x1);
else
((struct svm_parameter *) xpar_parameters->options)->nu = xpar_parameters->x1;
((struct svm_parameter *) xpar_parameters->options)->gamma = pow(10, xpar_parameters->x2);
}
static struct xpar_crit_parameters *init_xpar(struct dataset *dset, FILE *fdbg)
{
struct xpar_crit_parameters *xpar_params;
xpar_params = calloc(1, sizeof(struct xpar_crit_parameters));
if (xpar_params)
{
xpar_params->options = calloc(1, sizeof(struct svm_parameter));
if (xpar_params->options)
{
xpar_params->dset = dset;
xpar_params->outdev = stdout;
xpar_params->fdbg = fdbg;
xpar_params->vid = 0;
xpar_params->eval = FLT_MAX;
xpar_params->method = PALG_SVM;
xpar_params->nmodels = 1;
svm_param_defaults(xpar_params->options);
}
else
vx_free(xpar_params);
}
return xpar_params;
}
/*
Calculate criterion function values in vertices of `smx'.
*/
static int calc_simplex(int type, float *fval, float **smx, int ndim, void *xpar_parameters,
int triangle, float *ftrv, int *errc)
{
int i;
int status;
status = 0;
if (triangle == 0)
{
/*
First attempt, do the `base' triangle.
*/
if (type == C_SVC)
{
smx[0][0] = PCP_C_LEFT; /* C, log10-scale */
smx[1][0] = PCP_C_RIGHT;
smx[2][0] = PCP_C_MIDDLE;
}
else
{
smx[0][0] = PCP_NU_LEFT; /* nu */
smx[1][0] = get_svm_nu_max(tds->nd, tds->c);
smx[2][0] = 0.5*(smx[0][0]+smx[1][0]);
}
smx[0][1] = PCP_GAMMA_LEFT; /* gamma, log10-scale */
smx[1][1] = PCP_GAMMA_LEFT;
smx[2][1] = PCP_GAMMA_MIDDLE;
fval[0] = xpar_func(smx[0], ndim, 0, xpar_parameters, errc);
if (*errc == 0)
{
ftrv[0] = fval[0];
for (i = 1; (i < ndim+1) && !status; i++)
{
((struct xpar_crit_parameters *) xpar_parameters)->vid = i;
inverse_video();
fval[i] = xpar_func(smx[i], ndim, 0, xpar_parameters, errc);
if (*errc != 0)
status = -1;
else
ftrv[i] = fval[i];
}
}
else
status = -1;
}
else if (triangle == 1)
{
/*
Next. Simplex vertex 0 moves diagonally across, the other two
remain.
*/
if (type == C_SVC)
smx[0][0] = PCP_C_RIGHT;
else
smx[0][0] = get_svm_nu_max(tds->nd, tds->c);
smx[0][1] = PCP_GAMMA_RIGHT;
((struct xpar_crit_parameters *) xpar_parameters)->vid = 0;
fval[0] = xpar_func(smx[0], ndim, 0, xpar_parameters, errc);
if (*errc == 0)
{
ftrv[0] = fval[0];
fval[1] = ftrv[1];
fval[2] = ftrv[2];
}
else
status = -1;
}
else if (triangle == 2)
{
/*
Next. Simplex vertex 1 moves diagonally across, the other two
remain.
*/
if (type == C_SVC)
smx[1][0] = PCP_C_LEFT;
else
smx[1][0] = PCP_NU_LEFT;
smx[1][1] = PCP_GAMMA_RIGHT;
((struct xpar_crit_parameters *) xpar_parameters)->vid = 1;
fval[1] = xpar_func(smx[1], ndim, 0, xpar_parameters, errc);
if (*errc == 0)
{
ftrv[1] = fval[1];
fval[0] = ftrv[0];
fval[2] = ftrv[2];
}
else
status = -1;
}
else if (triangle == 3)
{
/*
Last triangle. Simplex vertex 0 moves diagonally across, the
other two remain.
*/
if (type == C_SVC)
smx[0][0] = PCP_C_LEFT;
else
smx[0][0] = PCP_NU_LEFT;
smx[0][1] = PCP_GAMMA_LEFT;
((struct xpar_crit_parameters *) xpar_parameters)->vid = 0;
fval[0] = xpar_func(smx[0], ndim, 0, xpar_parameters, errc);
if (*errc == 0)
{
ftrv[0] = fval[0];
fval[1] = ftrv[1];
fval[2] = ftrv[2];
}
else
status = -1;
}
return status;
}
/*
Optimize SVM parameters using Simplex algorithm.
*/
void pcp_svm_simplex(int *errc, int dbg, char **xname)
{
int ndim;
int itmax;
int i;
int iter;
int status;
int ex;
int type;
int done;
int triangle;
int dflt;
int idx;
float ftol = 1e-6;
float val;
char *msg;
float *fval;
float **smx;
FILE *outdev;
struct xpar_crit_parameters *xpar_parameters;
struct svm_parameter *svm_parm;
float ftrv[5];
FILE *fdbg = (FILE *) 0;
iter = 0;
if (dbg > 0)
fdbg = fopen(PCP_DBG, "a");
else
fdbg = (FILE *) 0;
ndim = 2;
status = smplx_alloc(&fval, &smx, ndim);
if (!status)
{
outdev = stdout;
xpar_parameters = init_xpar(tds, fdbg);
idx = ivec_min(tds->nd, tds->c);
svm_parm = calloc(1, sizeof(struct svm_parameter));
status = input_svm(outdev, tds->c, tds->nv, (struct svm_problem *) 0, svm_parm,
(char **) 0, (int *) 0, idx, 0, PCP_SVM_K_NONE, (int *) 0);
xpar_parameters->options = svm_parm;
input_xpar(tds, xpar_parameters);
dflt = PCP_DFLT_NITER;
ex = 1;
msg = malloc((PCP_QLEN+1)*sizeof(char));
snprintf(msg, PCP_QLEN+1, PCP_UMSG_MAXIT, dflt);
itmax = input_integer(stdin, outdev, msg, PCP_QLEN, &dflt, &ex, (int *) 0);
free(msg);
type = ((struct svm_parameter *) xpar_parameters->options)->svm_type;
/*
Compute feature transformations/subsets here, reuse them
later.
*/
if (xpar_parameters->dr_method != PDR_NONE)
status = compute_dr(xpar_parameters, xpar_parameters->dr_method,
xpar_parameters->idr, xpar_parameters->fscrit, fdbg, errc);
done = 0;
/*
We try four initial simplexes (triangles) until we execute
more than one iteration. ftrv holds the function values in
triangle vertices.
*/
triangle = 0;
unlink(PCP_MSL);
while (!done && !status)
{
inverse_video();
status = calc_simplex(type, fval, smx, ndim, (void *) xpar_parameters,
triangle, ftrv, errc);
if (!status)
{
inverse_video();
val = simplex(smx, ndim, fval, xpar_func, (void *) xpar_parameters,
itmax, ftol, &iter, errc);
if (*errc == 0)
{
if ((iter > 1) || (triangle == 3))
{
done = 1;
i = fvec_argmin(fval, ndim+1);
printf("\n");
viprint_line(4, 1, "Optimization completed in %7d iterations; minimum error rate: %7.2f%%.",
iter, xpar_parameters->eval);
if (type == C_SVC)
viprint_line(4, 1, "Optimal parameters: C = %12.5g; gamma = %12.5g",
pow(10, xpar_parameters->x1), pow(10, xpar_parameters->x2));
else
viprint_line(4, 1, "Optimal parameters: nu = %12.5g; gamma = %12.5g",
xpar_parameters->x1, pow(10, xpar_parameters->x2));
if (!teds)
pwait();
}
else
/*
If we finished in one iteration, that means no
progress has been made. Try one more time, with
a different initial simplex (i.e., next
triangle).
*/
triangle++;
}
else
status = -1;
}
}
/*
We have the optimal parameters - now apply to the test set, if
defined.
*/
if (teds)
{
printf("\n");
viprint_line(2, 1, "Applying the optimal classifier to the test dataset...");
set_optimal(xpar_parameters);
xtest_optimal(xpar_parameters, errc, xname);
}
xpar_free(xpar_parameters);
mx_free((void **) smx, ndim+1);
vx_free(fval);
}
}
static void pcp_grid_loop(int svm_type, int kernel_type, int *iter, float *grid_point,
float step, float gleft, float gstep, int nstep1, int nstep2,
int ndim, struct xpar_crit_parameters *xpar_parameters, int *errc)
{
int i;
int j;
float fval;
float pm; /* probability of error of the majority classifier */
float p_error;
float p_error_opt;
pm = (float) ivec_max(xpar_parameters->dset->nd, xpar_parameters->dset->c);
pm = 100.0*(1.0-pm/xpar_parameters->dset->nv);
for (i = 0; (i <= nstep1) && !*errc; i++) /* C/nu loop */
{
grid_point[0] += step;
grid_point[1] = gleft-gstep;
for (j = 0; (j <= nstep2) && !*errc; j++) /* other parameter loop */
{
grid_point[1] += gstep;
fval = xpar_func(grid_point, ndim, *iter, xpar_parameters, errc);
p_error = fval/pm;
p_error_opt = xpar_parameters->eval/pm;
if (!*errc)
{
if (kernel_type == PCP_SVM_K_LINEAR)
{
viprint_line(5, 1, "Current: C = %10.3g; error rate: %5.1f%%; p-error: %7.3f.",
pow(10, grid_point[0]), fval, p_error);
viprint_line(6, 1, "Optimal: C = %10.3g; error rate: %5.1f%%; p-error: %7.3f.",
pow(10, xpar_parameters->x1), xpar_parameters->eval, p_error_opt);
}
else if (kernel_type == PCP_SVM_K_RBF)
{
if (svm_type == C_SVC)
{
viprint_line(6, 1, "Current: C = %10.3g, gamma = %10.3g; error: %5.1f%%; p-error: %7.3f.",
pow(10, grid_point[0]), pow(10, grid_point[1]), fval, p_error);
viprint_line(6, 1, "Optimal: C = %10.3g, gamma = %10.3g; error: %5.1f%%; p-error: %7.3f.",
pow(10, xpar_parameters->x1), pow(10, xpar_parameters->x2),
xpar_parameters->eval, p_error_opt);
}
else
{
viprint_line(6, 1, "Current: nu = %10.3g, gamma = %10.3g; error: %5.1f%%; p-error: %7.3f.",
grid_point[0], pow(10, grid_point[1]), fval, p_error);
viprint_line(6, 1, "Optimal: nu = %10.3g, gamma = %10.3g; error: %5.1f%%; p-error: %7.3f.",
xpar_parameters->x1, pow(10, xpar_parameters->x2), xpar_parameters->eval,
fval, p_error_opt);
}
}
printf("\n");
(*iter)++;
}
}
}
}
/*
SVM model selection using grid (i.e., exhaustive) search.
*/
void pcp_svm_grid(int *errc, int dbg, char **xname)
{
int iter;
int ndim;
int svm_type;
int status;
int nstep1;
int nstep2;
int kernel_type;
float left;
float right;
float step;
float gleft;
float gright;
float gstep;
float numax;
float grid_point[2];
char *msg;
struct xpar_crit_parameters *xpar_parameters;
FILE *fdbg;
if (dbg > 0)
fdbg = fopen(PCP_DBG, "a");
else
fdbg = (FILE *) 0;
status = 0;
clear_screen();
cursor_on();
ndim = 2;
xpar_parameters = init_xpar(tds, fdbg);
msg = malloc((PCP_QLEN+1)*sizeof(char));
svm_type = input_svm_type(msg, stdout);
kernel_type = input_kernel(msg, stdout);
free(msg);
((struct svm_parameter *) xpar_parameters->options)->svm_type = svm_type;
((struct svm_parameter *) xpar_parameters->options)->kernel_type =
convert_kernel_type(kernel_type);
nstep1 = PCP_NSTEP;
if (kernel_type == PCP_SVM_K_RBF)
nstep2 = PCP_NSTEP;
else if (kernel_type == PCP_SVM_K_LINEAR)
{
ndim = 1;
nstep2 = 0;
}
if (svm_type == C_SVC)
input_class_costs(xpar_parameters->options, tds->c);
input_xpar(tds, xpar_parameters);
if (svm_type == C_SVC)
{
left = PCP_C_LEFT;
right = PCP_C_RIGHT;
}
else
{
left = PCP_NU_LEFT;
/*
Shrink the feasible interval for nu slightly, to avoid
violating the constraint due to round-off error in the nu loop
(the outer loop below, indexed by `i').
*/
numax = PCP_NU_FACTOR*get_svm_nu_xmax(tds->nd, tds->c, xpar_parameters->nxval);
right = numax;
}
step = (right-left)/PCP_NSTEP;
gleft = PCP_GAMMA_LEFT;
gright = PCP_GAMMA_RIGHT;
gstep = (gright-gleft)/PCP_NSTEP;
cursor_off();
inverse_video();
grid_point[0] = left-step;
*errc = 0;
unlink(PCP_MSL);
iter = 1;
/*
Compute optimal feature transformations/subsets once, reuse them
later.
*/
if (xpar_parameters->dr_method != PDR_NONE)
status = compute_dr(xpar_parameters, xpar_parameters->dr_method,
xpar_parameters->idr, xpar_parameters->fscrit, fdbg, errc);
if (!status)
{
pcp_grid_loop(svm_type, kernel_type, &iter, grid_point, step, gleft, gstep, nstep1,
nstep2, ndim, xpar_parameters, errc);
/*
Narrow down the search around the best pair of parameters
found so far. TBD: this code is tested and fully
functional. However it is presently disabled since it doesn't
seem to improve the predictions.
*/
if (!*errc && 0)
{
if (svm_type == NU_SVC)
{
nstep1 = PCP_RSTEP;
nstep2 = PCP_RSTEP;
left = (1-PCP_RFACTOR)*xpar_parameters->x1;
right = (1+PCP_RFACTOR)*xpar_parameters->x1;
step = (right-left)/nstep1;
if (right > numax)
right = numax;
}
else
{
/*
C-SVM is relatively insensitive to the choice of the C
parameter. Consequently, we don't tune it any
more. Instead, we increase the precision of tuning for
gamma (see `nstep2' variable).
*/
left = xpar_parameters->x1;
step = 0.0;
nstep1 = 0;
nstep2 = 5*PCP_RSTEP;
}
grid_point[0] = left-step;
gleft = (1-PCP_RFACTOR)*xpar_parameters->x2;
gright = (1+PCP_RFACTOR)*xpar_parameters->x2;
gstep = (gright-gleft)/nstep2;
pcp_grid_loop(svm_type, kernel_type, &iter, grid_point, step, gleft,
gstep, nstep1, nstep2, ndim, xpar_parameters, errc);
}
/*
We have the optimal parameters - now apply to the test set, if
defined.
*/
if (!*errc)
{
if (teds)
{
viprint_line(2, 1, "Applying the optimal classifier to the test dataset...");
set_optimal(xpar_parameters);
xtest_optimal(xpar_parameters, errc, xname);
}
else
pwait();
}
}
xpar_free(xpar_parameters);
}
syntax highlighted by Code2HTML, v. 0.9.1