diff --git a/.gitignore b/.gitignore index d4fd522..4a0091f 100644 --- a/.gitignore +++ b/.gitignore @@ -4,3 +4,4 @@ test*.jl deps/build.log +Manifest.toml diff --git a/.travis.yml b/.travis.yml index 5ebd08c..9834486 100644 --- a/.travis.yml +++ b/.travis.yml @@ -1,7 +1,7 @@ language: julia julia: - - "1.0" - "1.3" + - "1.4" - nightly os: - linux diff --git a/Project.toml b/Project.toml index 792e11b..3009590 100644 --- a/Project.toml +++ b/Project.toml @@ -1,19 +1,17 @@ name = "LIBSVM" uuid = "b1bec4e5-fd48-53fe-b0cb-9723c09d164b" -version = "0.4.0" +version = "0.5.0" [deps] -Compat = "34da2185-b29b-5c13-b0c7-acf172513d20" LIBLINEAR = "2d691ee1-e668-5016-a719-b2531b85e0f5" -Libdl = "8f399da3-3557-5675-b5ff-fb832c97cbdb" ScikitLearnBase = "6e75b9c4-186b-50bd-896f-2d2496a4843e" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" +libsvm_jll = "08558c22-525a-5d2a-acf6-0ac6658ffce4" [compat] -Compat = "2, 3" LIBLINEAR = "0.5" ScikitLearnBase = "0.5" -julia = "1" +julia = "1.3" [extras] DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab" diff --git a/appveyor.yml b/appveyor.yml index f01e29d..f56454f 100644 --- a/appveyor.yml +++ b/appveyor.yml @@ -1,6 +1,6 @@ environment: matrix: - - julia_version: "1.0" + - julia_version: "1.4" - julia_version: "1.3" # - julia_version: nightly diff --git a/deps/build.jl b/deps/build.jl deleted file mode 100644 index 9dbb0c2..0000000 --- a/deps/build.jl +++ /dev/null @@ -1,14 +0,0 @@ - -if Sys.iswindows() - lib = joinpath(joinpath(dirname(@__FILE__), "libsvm.dll")) - @info("Downloading LIBSVM binary") - if Sys.WORD_SIZE == 64 - download("http://web.ics.purdue.edu/~finej/libsvm-3.22_1.dll", lib) - else - download("http://web.ics.purdue.edu/~finej/libsvm32-3.22_1.dll", lib) - end -else - cd(joinpath(dirname(@__FILE__), "libsvm-3.22")) - run(`make lib`) - run(`mv libsvm.so.2 ../libsvm.so.2`) -end diff --git a/deps/libsvm-3.22/COPYRIGHT b/deps/libsvm-3.22/COPYRIGHT deleted file mode 100644 index 5fe2f22..0000000 --- a/deps/libsvm-3.22/COPYRIGHT +++ /dev/null @@ -1,31 +0,0 @@ - -Copyright (c) 2000-2014 Chih-Chung Chang and Chih-Jen Lin -All rights reserved. - -Redistribution and use in source and binary forms, with or without -modification, are permitted provided that the following conditions -are met: - -1. Redistributions of source code must retain the above copyright -notice, this list of conditions and the following disclaimer. - -2. Redistributions in binary form must reproduce the above copyright -notice, this list of conditions and the following disclaimer in the -documentation and/or other materials provided with the distribution. - -3. Neither name of copyright holders nor the names of its contributors -may be used to endorse or promote products derived from this software -without specific prior written permission. - - -THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS -``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT -LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR -A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR -CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, -EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, -PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR -PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF -LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING -NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS -SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. diff --git a/deps/libsvm-3.22/Makefile b/deps/libsvm-3.22/Makefile deleted file mode 100644 index ddf2aab..0000000 --- a/deps/libsvm-3.22/Makefile +++ /dev/null @@ -1,47 +0,0 @@ -CXX ?= g++ -CFLAGS_normal = -Wall -Wconversion -O3 -fPIC -DENABLEOPENMP -fopenmp -CFLAGS_fallback = -Wall -Wconversion -O3 -DDISABLEOPENMP -fPIC -SHVER = 2 -OS = $(shell uname) - -all: svm-train svm-predict svm-scale - -lib: - $(MAKE) lib_normal || $(MAKE) lib_fallback - -lib_normal: svm.o - if [ "$(OS)" = "Darwin" ]; then \ - SHARED_LIB_FLAG="-dynamiclib -Wl,-install_name,libsvm.so.$(SHVER)"; \ - else \ - SHARED_LIB_FLAG="-shared -Wl,-soname,libsvm.so.$(SHVER)"; \ - fi; \ - $(CXX) $${SHARED_LIB_FLAG} -DENABLEOPENMP -fopenmp svm.o -o libsvm.so.$(SHVER) - -lib_fallback: svm.o - if [ "$(OS)" = "Darwin" ]; then \ - SHARED_LIB_FLAG="-dynamiclib -Wl,-install_name,libsvm.so.$(SHVER)"; \ - else \ - SHARED_LIB_FLAG="-shared -Wl,-soname,libsvm.so.$(SHVER)"; \ - fi; \ - $(CXX) $${SHARED_LIB_FLAG} -DDISABLEOPENMP svm.o -o libsvm.so.$(SHVER) - -svm-predict: svm-predict.c svm.o - $(CXX) $(CFLAGS) svm-predict.c svm.o -o svm-predict -lm - -svm-train: svm-train.c svm.o - $(CXX) $(CFLAGS) svm-train.c svm.o -o svm-train -lm - -svm-scale: svm-scale.c - $(CXX) $(CFLAGS) svm-scale.c -o svm-scale - -svm.o: - $(MAKE) svm.o_normal || $(MAKE) svm.o_fallback - -svm.o_normal: svm.cpp svm.h - $(CXX) $(CFLAGS_normal) -c svm.cpp - -svm.o_fallback: svm.cpp svm.h - $(CXX) $(CFLAGS_fallblack) -c svm.cpp - -clean: - rm -f *~ svm.o svm-train svm-predict svm-scale libsvm.so.$(SHVER) diff --git a/deps/libsvm-3.22/Makefile.win32 b/deps/libsvm-3.22/Makefile.win32 deleted file mode 100644 index fd775fa..0000000 --- a/deps/libsvm-3.22/Makefile.win32 +++ /dev/null @@ -1,10 +0,0 @@ -CXX ?= cl.exe -CFLAGS = nologo /O2 /EHsc /I. /D _CRT_SECURE_NO_DEPRECATE /openmp - -all: lib - -lib: svm.cpp svm.h svm.def - $(CXX) $(CFLAGS) -LD svm.cpp -Felibsvm -link -DEF:svm.def - -clean: - rm -f *~ svm.obj libsvm.dll libsvm.lib libsvm.exp diff --git a/deps/libsvm-3.22/Makefile.win64 b/deps/libsvm-3.22/Makefile.win64 deleted file mode 100644 index 3285808..0000000 --- a/deps/libsvm-3.22/Makefile.win64 +++ /dev/null @@ -1,10 +0,0 @@ -CXX ?= cl.exe -CFLAGS = nologo /O2 /EHsc /I. /D _WIN64 /D _CRT_SECURE_NO_DEPRECATE /openmp - -all: lib - -lib: svm.cpp svm.h svm.def - $(CXX) $(CFLAGS) -LD svm.cpp -Felibsvm -link -DEF:svm.def - -clean: - rm -f *~ svm.obj libsvm.dll libsvm.lib libsvm.exp diff --git a/deps/libsvm-3.22/svm.cpp b/deps/libsvm-3.22/svm.cpp deleted file mode 100644 index 5175436..0000000 --- a/deps/libsvm-3.22/svm.cpp +++ /dev/null @@ -1,3203 +0,0 @@ -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include "svm.h" - -#ifdef ENABLEOPENMP - #include "omp.h" -#else - #define omp_get_max_threads() 0 - #define omp_set_num_threads(num_threads) void -#endif - -int libsvm_version = LIBSVM_VERSION; -typedef float Qfloat; -typedef signed char schar; -#ifndef min -template static inline T min(T x,T y) { return (x static inline T max(T x,T y) { return (x>y)?x:y; } -#endif -template static inline void swap(T& x, T& y) { T t=x; x=y; y=t; } -template static inline void clone(T*& dst, S* src, int n) -{ - dst = new T[n]; - memcpy((void *)dst,(void *)src,sizeof(T)*n); -} -static inline double powi(double base, int times) -{ - double tmp = base, ret = 1.0; - - for(int t=times; t>0; t/=2) - { - if(t%2==1) ret*=tmp; - tmp = tmp * tmp; - } - return ret; -} -#define INF HUGE_VAL -#define TAU 1e-12 -#define Malloc(type,n) (type *)malloc((n)*sizeof(type)) - -static void print_string_stdout(const char *s) -{ - fputs(s,stdout); - fflush(stdout); -} -static void (*svm_print_string) (const char *) = &print_string_stdout; -#if 1 -static void info(const char *fmt,...) -{ - char buf[BUFSIZ]; - va_list ap; - va_start(ap,fmt); - vsprintf(buf,fmt,ap); - va_end(ap); - (*svm_print_string)(buf); -} -#else -static void info(const char *fmt,...) {} -#endif - -// -// Kernel Cache -// -// l is the number of total data items -// size is the cache size limit in bytes -// -class Cache -{ -public: - Cache(int l,long int size); - ~Cache(); - - // request data [0,len) - // return some position p where [p,len) need to be filled - // (p >= len if nothing needs to be filled) - int get_data(const int index, Qfloat **data, int len); - void swap_index(int i, int j); -private: - int l; - long int size; - struct head_t - { - head_t *prev, *next; // a circular list - Qfloat *data; - int len; // data[0,len) is cached in this entry - }; - - head_t *head; - head_t lru_head; - void lru_delete(head_t *h); - void lru_insert(head_t *h); -}; - -Cache::Cache(int l_,long int size_):l(l_),size(size_) -{ - head = (head_t *)calloc(l,sizeof(head_t)); // initialized to 0 - size /= sizeof(Qfloat); - size -= l * sizeof(head_t) / sizeof(Qfloat); - size = max(size, 2 * (long int) l); // cache must be large enough for two columns - lru_head.next = lru_head.prev = &lru_head; -} - -Cache::~Cache() -{ - for(head_t *h = lru_head.next; h != &lru_head; h=h->next) - free(h->data); - free(head); -} - -void Cache::lru_delete(head_t *h) -{ - // delete from current location - h->prev->next = h->next; - h->next->prev = h->prev; -} - -void Cache::lru_insert(head_t *h) -{ - // insert to last position - h->next = &lru_head; - h->prev = lru_head.prev; - h->prev->next = h; - h->next->prev = h; -} - -int Cache::get_data(const int index, Qfloat **data, int len) -{ - head_t *h = &head[index]; - if(h->len) lru_delete(h); - int more = len - h->len; - - if(more > 0) - { - // free old space - while(size < more) - { - head_t *old = lru_head.next; - lru_delete(old); - free(old->data); - size += old->len; - old->data = 0; - old->len = 0; - } - - // allocate new space - h->data = (Qfloat *)realloc(h->data,sizeof(Qfloat)*len); - size -= more; - swap(h->len,len); - } - - lru_insert(h); - *data = h->data; - return len; -} - -void Cache::swap_index(int i, int j) -{ - if(i==j) return; - - if(head[i].len) lru_delete(&head[i]); - if(head[j].len) lru_delete(&head[j]); - swap(head[i].data,head[j].data); - swap(head[i].len,head[j].len); - if(head[i].len) lru_insert(&head[i]); - if(head[j].len) lru_insert(&head[j]); - - if(i>j) swap(i,j); - for(head_t *h = lru_head.next; h!=&lru_head; h=h->next) - { - if(h->len > i) - { - if(h->len > j) - swap(h->data[i],h->data[j]); - else - { - // give up - lru_delete(h); - free(h->data); - size += h->len; - h->data = 0; - h->len = 0; - } - } - } -} - -// -// Kernel evaluation -// -// the static method k_function is for doing single kernel evaluation -// the constructor of Kernel prepares to calculate the l*l kernel matrix -// the member function get_Q is for getting one column from the Q Matrix -// -class QMatrix { -public: - virtual Qfloat *get_Q(int column, int len) const = 0; - virtual double *get_QD() const = 0; - virtual void swap_index(int i, int j) const = 0; - virtual ~QMatrix() {} -}; - -class Kernel: public QMatrix { -public: - Kernel(int l, svm_node * const * x, const svm_parameter& param); - virtual ~Kernel(); - - static double k_function(const svm_node *x, const svm_node *y, - const svm_parameter& param); - virtual Qfloat *get_Q(int column, int len) const = 0; - virtual double *get_QD() const = 0; - virtual void swap_index(int i, int j) const // no so const... - { - swap(x[i],x[j]); - if(x_square) swap(x_square[i],x_square[j]); - } -protected: - - double (Kernel::*kernel_function)(int i, int j) const; - -private: - const svm_node **x; - double *x_square; - - // svm_parameter - const int kernel_type; - const int degree; - const double gamma; - const double coef0; - - static double dot(const svm_node *px, const svm_node *py); - double kernel_linear(int i, int j) const - { - return dot(x[i],x[j]); - } - double kernel_poly(int i, int j) const - { - return powi(gamma*dot(x[i],x[j])+coef0,degree); - } - double kernel_rbf(int i, int j) const - { - return exp(-gamma*(x_square[i]+x_square[j]-2*dot(x[i],x[j]))); - } - double kernel_sigmoid(int i, int j) const - { - return tanh(gamma*dot(x[i],x[j])+coef0); - } - double kernel_precomputed(int i, int j) const - { - return x[i][(int)(x[j][0].value)].value; - } -}; - -Kernel::Kernel(int l, svm_node * const * x_, const svm_parameter& param) -:kernel_type(param.kernel_type), degree(param.degree), - gamma(param.gamma), coef0(param.coef0) -{ - switch(kernel_type) - { - case LINEAR: - kernel_function = &Kernel::kernel_linear; - break; - case POLY: - kernel_function = &Kernel::kernel_poly; - break; - case RBF: - kernel_function = &Kernel::kernel_rbf; - break; - case SIGMOID: - kernel_function = &Kernel::kernel_sigmoid; - break; - case PRECOMPUTED: - kernel_function = &Kernel::kernel_precomputed; - break; - } - - clone(x,x_,l); - - if(kernel_type == RBF) - { - x_square = new double[l]; - for(int i=0;iindex != -1 && py->index != -1) - { - if(px->index == py->index) - { - sum += px->value * py->value; - ++px; - ++py; - } - else - { - if(px->index > py->index) - ++py; - else - ++px; - } - } - return sum; -} - -double Kernel::k_function(const svm_node *x, const svm_node *y, - const svm_parameter& param) -{ - switch(param.kernel_type) - { - case LINEAR: - return dot(x,y); - case POLY: - return powi(param.gamma*dot(x,y)+param.coef0,param.degree); - case RBF: - { - double sum = 0; - while(x->index != -1 && y->index !=-1) - { - if(x->index == y->index) - { - double d = x->value - y->value; - sum += d*d; - ++x; - ++y; - } - else - { - if(x->index > y->index) - { - sum += y->value * y->value; - ++y; - } - else - { - sum += x->value * x->value; - ++x; - } - } - } - - while(x->index != -1) - { - sum += x->value * x->value; - ++x; - } - - while(y->index != -1) - { - sum += y->value * y->value; - ++y; - } - - return exp(-param.gamma*sum); - } - case SIGMOID: - return tanh(param.gamma*dot(x,y)+param.coef0); - case PRECOMPUTED: //x: test (validation), y: SV - return x[(int)(y->value)].value; - default: - return 0; // Unreachable - } -} - -// An SMO algorithm in Fan et al., JMLR 6(2005), p. 1889--1918 -// Solves: -// -// min 0.5(\alpha^T Q \alpha) + p^T \alpha -// -// y^T \alpha = \delta -// y_i = +1 or -1 -// 0 <= alpha_i <= Cp for y_i = 1 -// 0 <= alpha_i <= Cn for y_i = -1 -// -// Given: -// -// Q, p, y, Cp, Cn, and an initial feasible point \alpha -// l is the size of vectors and matrices -// eps is the stopping tolerance -// -// solution will be put in \alpha, objective value will be put in obj -// -class Solver { -public: - Solver() {}; - virtual ~Solver() {}; - - struct SolutionInfo { - double obj; - double rho; - double upper_bound_p; - double upper_bound_n; - double r; // for Solver_NU - }; - - void Solve(int l, const QMatrix& Q, const double *p_, const schar *y_, - double *alpha_, double Cp, double Cn, double eps, - SolutionInfo* si, int shrinking); -protected: - int active_size; - schar *y; - double *G; // gradient of objective function - enum { LOWER_BOUND, UPPER_BOUND, FREE }; - char *alpha_status; // LOWER_BOUND, UPPER_BOUND, FREE - double *alpha; - const QMatrix *Q; - const double *QD; - double eps; - double Cp,Cn; - double *p; - int *active_set; - double *G_bar; // gradient, if we treat free variables as 0 - int l; - bool unshrink; // XXX - - double get_C(int i) - { - return (y[i] > 0)? Cp : Cn; - } - void update_alpha_status(int i) - { - if(alpha[i] >= get_C(i)) - alpha_status[i] = UPPER_BOUND; - else if(alpha[i] <= 0) - alpha_status[i] = LOWER_BOUND; - else alpha_status[i] = FREE; - } - bool is_upper_bound(int i) { return alpha_status[i] == UPPER_BOUND; } - bool is_lower_bound(int i) { return alpha_status[i] == LOWER_BOUND; } - bool is_free(int i) { return alpha_status[i] == FREE; } - void swap_index(int i, int j); - void reconstruct_gradient(); - virtual int select_working_set(int &i, int &j); - virtual double calculate_rho(); - virtual void do_shrinking(); -private: - bool be_shrunk(int i, double Gmax1, double Gmax2); -}; - -void Solver::swap_index(int i, int j) -{ - Q->swap_index(i,j); - swap(y[i],y[j]); - swap(G[i],G[j]); - swap(alpha_status[i],alpha_status[j]); - swap(alpha[i],alpha[j]); - swap(p[i],p[j]); - swap(active_set[i],active_set[j]); - swap(G_bar[i],G_bar[j]); -} - -void Solver::reconstruct_gradient() -{ - // reconstruct inactive elements of G from G_bar and free variables - - if(active_size == l) return; - - int i,j; - int nr_free = 0; - - for(j=active_size;j 2*active_size*(l-active_size)) - { - for(i=active_size;iget_Q(i,active_size); - for(j=0;jget_Q(i,l); - double alpha_i = alpha[i]; - for(j=active_size;jl = l; - this->Q = &Q; - QD=Q.get_QD(); - clone(p, p_,l); - clone(y, y_,l); - clone(alpha,alpha_,l); - this->Cp = Cp; - this->Cn = Cn; - this->eps = eps; - unshrink = false; - - // initialize alpha_status - { - alpha_status = new char[l]; - for(int i=0;iINT_MAX/100 ? INT_MAX : 100*l); - int counter = min(l,1000)+1; - - while(iter < max_iter) - { - // show progress and do shrinking - - if(--counter == 0) - { - counter = min(l,1000); - if(shrinking) do_shrinking(); - info("."); - } - - int i,j; - if(select_working_set(i,j)!=0) - { - // reconstruct the whole gradient - reconstruct_gradient(); - // reset active set size and check - active_size = l; - info("*"); - if(select_working_set(i,j)!=0) - break; - else - counter = 1; // do shrinking next iteration - } - - ++iter; - - // update alpha[i] and alpha[j], handle bounds carefully - - const Qfloat *Q_i = Q.get_Q(i,active_size); - const Qfloat *Q_j = Q.get_Q(j,active_size); - - double C_i = get_C(i); - double C_j = get_C(j); - - double old_alpha_i = alpha[i]; - double old_alpha_j = alpha[j]; - - if(y[i]!=y[j]) - { - double quad_coef = QD[i]+QD[j]+2*Q_i[j]; - if (quad_coef <= 0) - quad_coef = TAU; - double delta = (-G[i]-G[j])/quad_coef; - double diff = alpha[i] - alpha[j]; - alpha[i] += delta; - alpha[j] += delta; - - if(diff > 0) - { - if(alpha[j] < 0) - { - alpha[j] = 0; - alpha[i] = diff; - } - } - else - { - if(alpha[i] < 0) - { - alpha[i] = 0; - alpha[j] = -diff; - } - } - if(diff > C_i - C_j) - { - if(alpha[i] > C_i) - { - alpha[i] = C_i; - alpha[j] = C_i - diff; - } - } - else - { - if(alpha[j] > C_j) - { - alpha[j] = C_j; - alpha[i] = C_j + diff; - } - } - } - else - { - double quad_coef = QD[i]+QD[j]-2*Q_i[j]; - if (quad_coef <= 0) - quad_coef = TAU; - double delta = (G[i]-G[j])/quad_coef; - double sum = alpha[i] + alpha[j]; - alpha[i] -= delta; - alpha[j] += delta; - - if(sum > C_i) - { - if(alpha[i] > C_i) - { - alpha[i] = C_i; - alpha[j] = sum - C_i; - } - } - else - { - if(alpha[j] < 0) - { - alpha[j] = 0; - alpha[i] = sum; - } - } - if(sum > C_j) - { - if(alpha[j] > C_j) - { - alpha[j] = C_j; - alpha[i] = sum - C_j; - } - } - else - { - if(alpha[i] < 0) - { - alpha[i] = 0; - alpha[j] = sum; - } - } - } - - // update G - - double delta_alpha_i = alpha[i] - old_alpha_i; - double delta_alpha_j = alpha[j] - old_alpha_j; - - for(int k=0;k= max_iter) - { - if(active_size < l) - { - // reconstruct the whole gradient to calculate objective value - reconstruct_gradient(); - active_size = l; - info("*"); - } - fprintf(stderr,"\nWARNING: reaching max number of iterations\n"); - } - - // calculate rho - - si->rho = calculate_rho(); - - // calculate objective value - { - double v = 0; - int i; - for(i=0;iobj = v/2; - } - - // put back the solution - { - for(int i=0;iupper_bound_p = Cp; - si->upper_bound_n = Cn; - - info("\noptimization finished, #iter = %d\n",iter); - - delete[] p; - delete[] y; - delete[] alpha; - delete[] alpha_status; - delete[] active_set; - delete[] G; - delete[] G_bar; -} - -// return 1 if already optimal, return 0 otherwise -int Solver::select_working_set(int &out_i, int &out_j) -{ - // return i,j such that - // i: maximizes -y_i * grad(f)_i, i in I_up(\alpha) - // j: minimizes the decrease of obj value - // (if quadratic coefficeint <= 0, replace it with tau) - // -y_j*grad(f)_j < -y_i*grad(f)_i, j in I_low(\alpha) - - double Gmax = -INF; - double Gmax2 = -INF; - int Gmax_idx = -1; - int Gmin_idx = -1; - double obj_diff_min = INF; - - for(int t=0;t= Gmax) - { - Gmax = -G[t]; - Gmax_idx = t; - } - } - else - { - if(!is_lower_bound(t)) - if(G[t] >= Gmax) - { - Gmax = G[t]; - Gmax_idx = t; - } - } - - int i = Gmax_idx; - const Qfloat *Q_i = NULL; - if(i != -1) // NULL Q_i not accessed: Gmax=-INF if i=-1 - Q_i = Q->get_Q(i,active_size); - - for(int j=0;j= Gmax2) - Gmax2 = G[j]; - if (grad_diff > 0) - { - double obj_diff; - double quad_coef = QD[i]+QD[j]-2.0*y[i]*Q_i[j]; - if (quad_coef > 0) - obj_diff = -(grad_diff*grad_diff)/quad_coef; - else - obj_diff = -(grad_diff*grad_diff)/TAU; - - if (obj_diff <= obj_diff_min) - { - Gmin_idx=j; - obj_diff_min = obj_diff; - } - } - } - } - else - { - if (!is_upper_bound(j)) - { - double grad_diff= Gmax-G[j]; - if (-G[j] >= Gmax2) - Gmax2 = -G[j]; - if (grad_diff > 0) - { - double obj_diff; - double quad_coef = QD[i]+QD[j]+2.0*y[i]*Q_i[j]; - if (quad_coef > 0) - obj_diff = -(grad_diff*grad_diff)/quad_coef; - else - obj_diff = -(grad_diff*grad_diff)/TAU; - - if (obj_diff <= obj_diff_min) - { - Gmin_idx=j; - obj_diff_min = obj_diff; - } - } - } - } - } - - if(Gmax+Gmax2 < eps || Gmin_idx == -1) - return 1; - - out_i = Gmax_idx; - out_j = Gmin_idx; - return 0; -} - -bool Solver::be_shrunk(int i, double Gmax1, double Gmax2) -{ - if(is_upper_bound(i)) - { - if(y[i]==+1) - return(-G[i] > Gmax1); - else - return(-G[i] > Gmax2); - } - else if(is_lower_bound(i)) - { - if(y[i]==+1) - return(G[i] > Gmax2); - else - return(G[i] > Gmax1); - } - else - return(false); -} - -void Solver::do_shrinking() -{ - int i; - double Gmax1 = -INF; // max { -y_i * grad(f)_i | i in I_up(\alpha) } - double Gmax2 = -INF; // max { y_i * grad(f)_i | i in I_low(\alpha) } - - // find maximal violating pair first - for(i=0;i= Gmax1) - Gmax1 = -G[i]; - } - if(!is_lower_bound(i)) - { - if(G[i] >= Gmax2) - Gmax2 = G[i]; - } - } - else - { - if(!is_upper_bound(i)) - { - if(-G[i] >= Gmax2) - Gmax2 = -G[i]; - } - if(!is_lower_bound(i)) - { - if(G[i] >= Gmax1) - Gmax1 = G[i]; - } - } - } - - if(unshrink == false && Gmax1 + Gmax2 <= eps*10) - { - unshrink = true; - reconstruct_gradient(); - active_size = l; - info("*"); - } - - for(i=0;i i) - { - if (!be_shrunk(active_size, Gmax1, Gmax2)) - { - swap_index(i,active_size); - break; - } - active_size--; - } - } -} - -double Solver::calculate_rho() -{ - double r; - int nr_free = 0; - double ub = INF, lb = -INF, sum_free = 0; - for(int i=0;i0) - r = sum_free/nr_free; - else - r = (ub+lb)/2; - - return r; -} - -// -// Solver for nu-svm classification and regression -// -// additional constraint: e^T \alpha = constant -// -class Solver_NU: public Solver -{ -public: - Solver_NU() {} - void Solve(int l, const QMatrix& Q, const double *p, const schar *y, - double *alpha, double Cp, double Cn, double eps, - SolutionInfo* si, int shrinking) - { - this->si = si; - Solver::Solve(l,Q,p,y,alpha,Cp,Cn,eps,si,shrinking); - } -private: - SolutionInfo *si; - int select_working_set(int &i, int &j); - double calculate_rho(); - bool be_shrunk(int i, double Gmax1, double Gmax2, double Gmax3, double Gmax4); - void do_shrinking(); -}; - -// return 1 if already optimal, return 0 otherwise -int Solver_NU::select_working_set(int &out_i, int &out_j) -{ - // return i,j such that y_i = y_j and - // i: maximizes -y_i * grad(f)_i, i in I_up(\alpha) - // j: minimizes the decrease of obj value - // (if quadratic coefficeint <= 0, replace it with tau) - // -y_j*grad(f)_j < -y_i*grad(f)_i, j in I_low(\alpha) - - double Gmaxp = -INF; - double Gmaxp2 = -INF; - int Gmaxp_idx = -1; - - double Gmaxn = -INF; - double Gmaxn2 = -INF; - int Gmaxn_idx = -1; - - int Gmin_idx = -1; - double obj_diff_min = INF; - - for(int t=0;t= Gmaxp) - { - Gmaxp = -G[t]; - Gmaxp_idx = t; - } - } - else - { - if(!is_lower_bound(t)) - if(G[t] >= Gmaxn) - { - Gmaxn = G[t]; - Gmaxn_idx = t; - } - } - - int ip = Gmaxp_idx; - int in = Gmaxn_idx; - const Qfloat *Q_ip = NULL; - const Qfloat *Q_in = NULL; - if(ip != -1) // NULL Q_ip not accessed: Gmaxp=-INF if ip=-1 - Q_ip = Q->get_Q(ip,active_size); - if(in != -1) - Q_in = Q->get_Q(in,active_size); - - for(int j=0;j= Gmaxp2) - Gmaxp2 = G[j]; - if (grad_diff > 0) - { - double obj_diff; - double quad_coef = QD[ip]+QD[j]-2*Q_ip[j]; - if (quad_coef > 0) - obj_diff = -(grad_diff*grad_diff)/quad_coef; - else - obj_diff = -(grad_diff*grad_diff)/TAU; - - if (obj_diff <= obj_diff_min) - { - Gmin_idx=j; - obj_diff_min = obj_diff; - } - } - } - } - else - { - if (!is_upper_bound(j)) - { - double grad_diff=Gmaxn-G[j]; - if (-G[j] >= Gmaxn2) - Gmaxn2 = -G[j]; - if (grad_diff > 0) - { - double obj_diff; - double quad_coef = QD[in]+QD[j]-2*Q_in[j]; - if (quad_coef > 0) - obj_diff = -(grad_diff*grad_diff)/quad_coef; - else - obj_diff = -(grad_diff*grad_diff)/TAU; - - if (obj_diff <= obj_diff_min) - { - Gmin_idx=j; - obj_diff_min = obj_diff; - } - } - } - } - } - - if(max(Gmaxp+Gmaxp2,Gmaxn+Gmaxn2) < eps || Gmin_idx == -1) - return 1; - - if (y[Gmin_idx] == +1) - out_i = Gmaxp_idx; - else - out_i = Gmaxn_idx; - out_j = Gmin_idx; - - return 0; -} - -bool Solver_NU::be_shrunk(int i, double Gmax1, double Gmax2, double Gmax3, double Gmax4) -{ - if(is_upper_bound(i)) - { - if(y[i]==+1) - return(-G[i] > Gmax1); - else - return(-G[i] > Gmax4); - } - else if(is_lower_bound(i)) - { - if(y[i]==+1) - return(G[i] > Gmax2); - else - return(G[i] > Gmax3); - } - else - return(false); -} - -void Solver_NU::do_shrinking() -{ - double Gmax1 = -INF; // max { -y_i * grad(f)_i | y_i = +1, i in I_up(\alpha) } - double Gmax2 = -INF; // max { y_i * grad(f)_i | y_i = +1, i in I_low(\alpha) } - double Gmax3 = -INF; // max { -y_i * grad(f)_i | y_i = -1, i in I_up(\alpha) } - double Gmax4 = -INF; // max { y_i * grad(f)_i | y_i = -1, i in I_low(\alpha) } - - // find maximal violating pair first - int i; - for(i=0;i Gmax1) Gmax1 = -G[i]; - } - else if(-G[i] > Gmax4) Gmax4 = -G[i]; - } - if(!is_lower_bound(i)) - { - if(y[i]==+1) - { - if(G[i] > Gmax2) Gmax2 = G[i]; - } - else if(G[i] > Gmax3) Gmax3 = G[i]; - } - } - - if(unshrink == false && max(Gmax1+Gmax2,Gmax3+Gmax4) <= eps*10) - { - unshrink = true; - reconstruct_gradient(); - active_size = l; - } - - for(i=0;i i) - { - if (!be_shrunk(active_size, Gmax1, Gmax2, Gmax3, Gmax4)) - { - swap_index(i,active_size); - break; - } - active_size--; - } - } -} - -double Solver_NU::calculate_rho() -{ - int nr_free1 = 0,nr_free2 = 0; - double ub1 = INF, ub2 = INF; - double lb1 = -INF, lb2 = -INF; - double sum_free1 = 0, sum_free2 = 0; - - for(int i=0;i 0) - r1 = sum_free1/nr_free1; - else - r1 = (ub1+lb1)/2; - - if(nr_free2 > 0) - r2 = sum_free2/nr_free2; - else - r2 = (ub2+lb2)/2; - - si->r = (r1+r2)/2; - return (r1-r2)/2; -} - -// -// Q matrices for various formulations -// -class SVC_Q: public Kernel -{ -public: - SVC_Q(const svm_problem& prob, const svm_parameter& param, const schar *y_) - :Kernel(prob.l, prob.x, param) - { - clone(y,y_,prob.l); - cache = new Cache(prob.l,(long int)(param.cache_size*(1<<20))); - QD = new double[prob.l]; - for(int i=0;i*kernel_function)(i,i); - } - - Qfloat *get_Q(int i, int len) const - { - Qfloat *data; - int start, j; - if((start = cache->get_data(i,&data,len)) < len) - { - #pragma omp parallel for private(j) schedule(guided) - for(j=start;j*kernel_function)(i,j)); - } - return data; - } - - double *get_QD() const - { - return QD; - } - - void swap_index(int i, int j) const - { - cache->swap_index(i,j); - Kernel::swap_index(i,j); - swap(y[i],y[j]); - swap(QD[i],QD[j]); - } - - ~SVC_Q() - { - delete[] y; - delete cache; - delete[] QD; - } -private: - schar *y; - Cache *cache; - double *QD; -}; - -class ONE_CLASS_Q: public Kernel -{ -public: - ONE_CLASS_Q(const svm_problem& prob, const svm_parameter& param) - :Kernel(prob.l, prob.x, param) - { - cache = new Cache(prob.l,(long int)(param.cache_size*(1<<20))); - QD = new double[prob.l]; - for(int i=0;i*kernel_function)(i,i); - } - - Qfloat *get_Q(int i, int len) const - { - Qfloat *data; - int start, j; - if((start = cache->get_data(i,&data,len)) < len) - { - for(j=start;j*kernel_function)(i,j); - } - return data; - } - - double *get_QD() const - { - return QD; - } - - void swap_index(int i, int j) const - { - cache->swap_index(i,j); - Kernel::swap_index(i,j); - swap(QD[i],QD[j]); - } - - ~ONE_CLASS_Q() - { - delete cache; - delete[] QD; - } -private: - Cache *cache; - double *QD; -}; - -class SVR_Q: public Kernel -{ -public: - SVR_Q(const svm_problem& prob, const svm_parameter& param) - :Kernel(prob.l, prob.x, param) - { - l = prob.l; - cache = new Cache(l,(long int)(param.cache_size*(1<<20))); - QD = new double[2*l]; - sign = new schar[2*l]; - index = new int[2*l]; - for(int k=0;k*kernel_function)(k,k); - QD[k+l] = QD[k]; - } - buffer[0] = new Qfloat[2*l]; - buffer[1] = new Qfloat[2*l]; - next_buffer = 0; - } - - void swap_index(int i, int j) const - { - swap(sign[i],sign[j]); - swap(index[i],index[j]); - swap(QD[i],QD[j]); - } - - Qfloat *get_Q(int i, int len) const - { - Qfloat *data; - int j, real_i = index[i]; - if(cache->get_data(real_i,&data,l) < l) - { - #pragma omp parallel for private(i) schedule(guided) - for(j=0;j*kernel_function)(real_i,j); - } - - // reorder and copy - Qfloat *buf = buffer[next_buffer]; - next_buffer = 1 - next_buffer; - schar si = sign[i]; - for(j=0;jl; - double *minus_ones = new double[l]; - schar *y = new schar[l]; - - int i; - - for(i=0;iy[i] > 0) y[i] = +1; else y[i] = -1; - } - - Solver s; - s.Solve(l, SVC_Q(*prob,*param,y), minus_ones, y, - alpha, Cp, Cn, param->eps, si, param->shrinking); - - double sum_alpha=0; - for(i=0;il)); - - for(i=0;il; - double nu = param->nu; - - schar *y = new schar[l]; - - for(i=0;iy[i]>0) - y[i] = +1; - else - y[i] = -1; - - double sum_pos = nu*l/2; - double sum_neg = nu*l/2; - - for(i=0;ieps, si, param->shrinking); - double r = si->r; - - info("C = %f\n",1/r); - - for(i=0;irho /= r; - si->obj /= (r*r); - si->upper_bound_p = 1/r; - si->upper_bound_n = 1/r; - - delete[] y; - delete[] zeros; -} - -static void solve_one_class( - const svm_problem *prob, const svm_parameter *param, - double *alpha, Solver::SolutionInfo* si) -{ - int l = prob->l; - double *zeros = new double[l]; - schar *ones = new schar[l]; - int i; - - int n = (int)(param->nu*prob->l); // # of alpha's at upper bound - - for(i=0;il) - alpha[n] = param->nu * prob->l - n; - for(i=n+1;ieps, si, param->shrinking); - - delete[] zeros; - delete[] ones; -} - -static void solve_epsilon_svr( - const svm_problem *prob, const svm_parameter *param, - double *alpha, Solver::SolutionInfo* si) -{ - int l = prob->l; - double *alpha2 = new double[2*l]; - double *linear_term = new double[2*l]; - schar *y = new schar[2*l]; - int i; - - for(i=0;ip - prob->y[i]; - y[i] = 1; - - alpha2[i+l] = 0; - linear_term[i+l] = param->p + prob->y[i]; - y[i+l] = -1; - } - - Solver s; - s.Solve(2*l, SVR_Q(*prob,*param), linear_term, y, - alpha2, param->C, param->C, param->eps, si, param->shrinking); - - double sum_alpha = 0; - for(i=0;iC*l)); - - delete[] alpha2; - delete[] linear_term; - delete[] y; -} - -static void solve_nu_svr( - const svm_problem *prob, const svm_parameter *param, - double *alpha, Solver::SolutionInfo* si) -{ - int l = prob->l; - double C = param->C; - double *alpha2 = new double[2*l]; - double *linear_term = new double[2*l]; - schar *y = new schar[2*l]; - int i; - - double sum = C * param->nu * l / 2; - for(i=0;iy[i]; - y[i] = 1; - - linear_term[i+l] = prob->y[i]; - y[i+l] = -1; - } - - Solver_NU s; - s.Solve(2*l, SVR_Q(*prob,*param), linear_term, y, - alpha2, C, C, param->eps, si, param->shrinking); - - info("epsilon = %f\n",-si->r); - - for(i=0;il); - Solver::SolutionInfo si; - switch(param->svm_type) - { - case C_SVC: - solve_c_svc(prob,param,alpha,&si,Cp,Cn); - break; - case NU_SVC: - solve_nu_svc(prob,param,alpha,&si); - break; - case ONE_CLASS: - solve_one_class(prob,param,alpha,&si); - break; - case EPSILON_SVR: - solve_epsilon_svr(prob,param,alpha,&si); - break; - case NU_SVR: - solve_nu_svr(prob,param,alpha,&si); - break; - } - - info("obj = %f, rho = %f\n",si.obj,si.rho); - - // output SVs - - int nSV = 0; - int nBSV = 0; - for(int i=0;il;i++) - { - if(fabs(alpha[i]) > 0) - { - ++nSV; - if(prob->y[i] > 0) - { - if(fabs(alpha[i]) >= si.upper_bound_p) - ++nBSV; - } - else - { - if(fabs(alpha[i]) >= si.upper_bound_n) - ++nBSV; - } - } - } - - info("nSV = %d, nBSV = %d\n",nSV,nBSV); - - decision_function f; - f.alpha = alpha; - f.rho = si.rho; - return f; -} - -// Platt's binary SVM Probablistic Output: an improvement from Lin et al. -static void sigmoid_train( - int l, const double *dec_values, const double *labels, - double& A, double& B) -{ - double prior1=0, prior0 = 0; - int i; - - for (i=0;i 0) prior1+=1; - else prior0+=1; - - int max_iter=100; // Maximal number of iterations - double min_step=1e-10; // Minimal step taken in line search - double sigma=1e-12; // For numerically strict PD of Hessian - double eps=1e-5; - double hiTarget=(prior1+1.0)/(prior1+2.0); - double loTarget=1/(prior0+2.0); - double *t=Malloc(double,l); - double fApB,p,q,h11,h22,h21,g1,g2,det,dA,dB,gd,stepsize; - double newA,newB,newf,d1,d2; - int iter; - - // Initial Point and Initial Fun Value - A=0.0; B=log((prior0+1.0)/(prior1+1.0)); - double fval = 0.0; - - for (i=0;i0) t[i]=hiTarget; - else t[i]=loTarget; - fApB = dec_values[i]*A+B; - if (fApB>=0) - fval += t[i]*fApB + log(1+exp(-fApB)); - else - fval += (t[i] - 1)*fApB +log(1+exp(fApB)); - } - for (iter=0;iter= 0) - { - p=exp(-fApB)/(1.0+exp(-fApB)); - q=1.0/(1.0+exp(-fApB)); - } - else - { - p=1.0/(1.0+exp(fApB)); - q=exp(fApB)/(1.0+exp(fApB)); - } - d2=p*q; - h11+=dec_values[i]*dec_values[i]*d2; - h22+=d2; - h21+=dec_values[i]*d2; - d1=t[i]-p; - g1+=dec_values[i]*d1; - g2+=d1; - } - - // Stopping Criteria - if (fabs(g1)= min_step) - { - newA = A + stepsize * dA; - newB = B + stepsize * dB; - - // New function value - newf = 0.0; - for (i=0;i= 0) - newf += t[i]*fApB + log(1+exp(-fApB)); - else - newf += (t[i] - 1)*fApB +log(1+exp(fApB)); - } - // Check sufficient decrease - if (newf=max_iter) - info("Reaching maximal iterations in two-class probability estimates\n"); - free(t); -} - -static double sigmoid_predict(double decision_value, double A, double B) -{ - double fApB = decision_value*A+B; - // 1-p used later; avoid catastrophic cancellation - if (fApB >= 0) - return exp(-fApB)/(1.0+exp(-fApB)); - else - return 1.0/(1+exp(fApB)) ; -} - -// Method 2 from the multiclass_prob paper by Wu, Lin, and Weng -static void multiclass_probability(int k, double **r, double *p) -{ - int t,j; - int iter = 0, max_iter=max(100,k); - double **Q=Malloc(double *,k); - double *Qp=Malloc(double,k); - double pQp, eps=0.005/k; - - for (t=0;tmax_error) - max_error=error; - } - if (max_error=max_iter) - info("Exceeds max_iter in multiclass_prob\n"); - for(t=0;tl); - double *dec_values = Malloc(double,prob->l); - - // random shuffle - for(i=0;il;i++) perm[i]=i; - for(i=0;il;i++) - { - int j = i+rand()%(prob->l-i); - swap(perm[i],perm[j]); - } - for(i=0;il/nr_fold; - int end = (i+1)*prob->l/nr_fold; - int j,k; - struct svm_problem subprob; - - subprob.l = prob->l-(end-begin); - subprob.x = Malloc(struct svm_node*,subprob.l); - subprob.y = Malloc(double,subprob.l); - - k=0; - for(j=0;jx[perm[j]]; - subprob.y[k] = prob->y[perm[j]]; - ++k; - } - for(j=end;jl;j++) - { - subprob.x[k] = prob->x[perm[j]]; - subprob.y[k] = prob->y[perm[j]]; - ++k; - } - int p_count=0,n_count=0; - for(j=0;j0) - p_count++; - else - n_count++; - - if(p_count==0 && n_count==0) - for(j=begin;j 0 && n_count == 0) - for(j=begin;j 0) - for(j=begin;jx[perm[j]],&(dec_values[perm[j]])); - // ensure +1 -1 order; reason not using CV subroutine - dec_values[perm[j]] *= submodel->label[0]; - } - svm_free_and_destroy_model(&submodel); - svm_destroy_param(&subparam); - } - free(subprob.x); - free(subprob.y); - } - sigmoid_train(prob->l,dec_values,prob->y,probA,probB); - free(dec_values); - free(perm); -} - -// Return parameter of a Laplace distribution -static double svm_svr_probability( - const svm_problem *prob, const svm_parameter *param) -{ - int i; - int nr_fold = 5; - double *ymv = Malloc(double,prob->l); - double mae = 0; - - svm_parameter newparam = *param; - newparam.probability = 0; - svm_cross_validation(prob,&newparam,nr_fold,ymv); - for(i=0;il;i++) - { - ymv[i]=prob->y[i]-ymv[i]; - mae += fabs(ymv[i]); - } - mae /= prob->l; - double std=sqrt(2*mae*mae); - int count=0; - mae=0; - for(i=0;il;i++) - if (fabs(ymv[i]) > 5*std) - count=count+1; - else - mae+=fabs(ymv[i]); - mae /= (prob->l-count); - info("Prob. model for test data: target value = predicted value + z,\nz: Laplace distribution e^(-|z|/sigma)/(2sigma),sigma= %g\n",mae); - free(ymv); - return mae; -} - - -// label: label name, start: begin of each class, count: #data of classes, perm: indices to the original data -// perm, length l, must be allocated before calling this subroutine -static void svm_group_classes(const svm_problem *prob, int *nr_class_ret, int **label_ret, int **start_ret, int **count_ret, int *perm) -{ - int l = prob->l; - int max_nr_class = 16; - int nr_class = 0; - int *label = Malloc(int,max_nr_class); - int *count = Malloc(int,max_nr_class); - int *data_label = Malloc(int,l); - int i; - - for(i=0;iy[i]; - int j; - for(j=0;jparam = *param; - model->free_sv = 0; // XXX - - if(param->svm_type == ONE_CLASS || - param->svm_type == EPSILON_SVR || - param->svm_type == NU_SVR) - { - // regression or one-class-svm - model->nr_class = 2; - model->label = NULL; - model->nSV = NULL; - model->probA = NULL; model->probB = NULL; - model->sv_coef = Malloc(double *,1); - - if(param->probability && - (param->svm_type == EPSILON_SVR || - param->svm_type == NU_SVR)) - { - model->probA = Malloc(double,1); - model->probA[0] = svm_svr_probability(prob,param); - } - - decision_function f = svm_train_one(prob,param,0,0); - model->rho = Malloc(double,1); - model->rho[0] = f.rho; - - int nSV = 0; - int i; - for(i=0;il;i++) - if(fabs(f.alpha[i]) > 0) ++nSV; - model->l = nSV; - model->SV = Malloc(svm_node *,nSV); - model->sv_coef[0] = Malloc(double,nSV); - model->sv_indices = Malloc(int,nSV); - int j = 0; - for(i=0;il;i++) - if(fabs(f.alpha[i]) > 0) - { - model->SV[j] = prob->x[i]; - model->sv_coef[0][j] = f.alpha[i]; - model->sv_indices[j] = i+1; - ++j; - } - - free(f.alpha); - } - else - { - // classification - int l = prob->l; - int nr_class; - int *label = NULL; - int *start = NULL; - int *count = NULL; - int *perm = Malloc(int,l); - - // group training data of the same class - svm_group_classes(prob,&nr_class,&label,&start,&count,perm); - if(nr_class == 1) - info("WARNING: training data in only one class. See README for details.\n"); - - svm_node **x = Malloc(svm_node *,l); - int i; - for(i=0;ix[perm[i]]; - - // calculate weighted C - - double *weighted_C = Malloc(double, nr_class); - for(i=0;iC; - for(i=0;inr_weight;i++) - { - int j; - for(j=0;jweight_label[i] == label[j]) - break; - if(j == nr_class) - fprintf(stderr,"WARNING: class label %d specified in weight is not found\n", param->weight_label[i]); - else - weighted_C[j] *= param->weight[i]; - } - - // train k*(k-1)/2 models - - bool *nonzero = Malloc(bool,l); - for(i=0;iprobability) - { - probA=Malloc(double,nr_class*(nr_class-1)/2); - probB=Malloc(double,nr_class*(nr_class-1)/2); - } - - int p = 0; - for(i=0;iprobability) - svm_binary_svc_probability(&sub_prob,param,weighted_C[i],weighted_C[j],probA[p],probB[p]); - - f[p] = svm_train_one(&sub_prob,param,weighted_C[i],weighted_C[j]); - for(k=0;k 0) - nonzero[si+k] = true; - for(k=0;k 0) - nonzero[sj+k] = true; - free(sub_prob.x); - free(sub_prob.y); - ++p; - } - - // build output - - model->nr_class = nr_class; - - model->label = Malloc(int,nr_class); - for(i=0;ilabel[i] = label[i]; - - model->rho = Malloc(double,nr_class*(nr_class-1)/2); - for(i=0;irho[i] = f[i].rho; - - if(param->probability) - { - model->probA = Malloc(double,nr_class*(nr_class-1)/2); - model->probB = Malloc(double,nr_class*(nr_class-1)/2); - for(i=0;iprobA[i] = probA[i]; - model->probB[i] = probB[i]; - } - } - else - { - model->probA=NULL; - model->probB=NULL; - } - - int total_sv = 0; - int *nz_count = Malloc(int,nr_class); - model->nSV = Malloc(int,nr_class); - for(i=0;inSV[i] = nSV; - nz_count[i] = nSV; - } - - info("Total nSV = %d\n",total_sv); - - model->l = total_sv; - model->SV = Malloc(svm_node *,total_sv); - model->sv_indices = Malloc(int,total_sv); - p = 0; - for(i=0;iSV[p] = x[i]; - model->sv_indices[p++] = perm[i] + 1; - } - - int *nz_start = Malloc(int,nr_class); - nz_start[0] = 0; - for(i=1;isv_coef = Malloc(double *,nr_class-1); - for(i=0;isv_coef[i] = Malloc(double,total_sv); - - p = 0; - for(i=0;isv_coef[j-1][q++] = f[p].alpha[k]; - q = nz_start[j]; - for(k=0;ksv_coef[i][q++] = f[p].alpha[ci+k]; - ++p; - } - - free(label); - free(probA); - free(probB); - free(count); - free(perm); - free(start); - free(x); - free(weighted_C); - free(nonzero); - for(i=0;il; - int *perm = Malloc(int,l); - int nr_class; - if (nr_fold > l) - { - nr_fold = l; - fprintf(stderr,"WARNING: # folds > # data. Will use # folds = # data instead (i.e., leave-one-out cross validation)\n"); - } - fold_start = Malloc(int,nr_fold+1); - // stratified cv may not give leave-one-out rate - // Each class to l folds -> some folds may have zero elements - if((param->svm_type == C_SVC || - param->svm_type == NU_SVC) && nr_fold < l) - { - int *start = NULL; - int *label = NULL; - int *count = NULL; - svm_group_classes(prob,&nr_class,&label,&start,&count,perm); - - // random shuffle and then data grouped by fold using the array perm - int *fold_count = Malloc(int,nr_fold); - int c; - int *index = Malloc(int,l); - for(i=0;ix[perm[j]]; - subprob.y[k] = prob->y[perm[j]]; - ++k; - } - for(j=end;jx[perm[j]]; - subprob.y[k] = prob->y[perm[j]]; - ++k; - } - struct svm_model *submodel = svm_train(&subprob,param); - if(param->probability && - (param->svm_type == C_SVC || param->svm_type == NU_SVC)) - { - double *prob_estimates=Malloc(double,svm_get_nr_class(submodel)); - for(j=begin;jx[perm[j]],prob_estimates); - free(prob_estimates); - } - else - for(j=begin;jx[perm[j]]); - svm_free_and_destroy_model(&submodel); - free(subprob.x); - free(subprob.y); - } - free(fold_start); - free(perm); -} - - -int svm_get_svm_type(const svm_model *model) -{ - return model->param.svm_type; -} - -int svm_get_nr_class(const svm_model *model) -{ - return model->nr_class; -} - -void svm_get_labels(const svm_model *model, int* label) -{ - if (model->label != NULL) - for(int i=0;inr_class;i++) - label[i] = model->label[i]; -} - -void svm_get_sv_indices(const svm_model *model, int* indices) -{ - if (model->sv_indices != NULL) - for(int i=0;il;i++) - indices[i] = model->sv_indices[i]; -} - -int svm_get_nr_sv(const svm_model *model) -{ - return model->l; -} - -double svm_get_svr_probability(const svm_model *model) -{ - if ((model->param.svm_type == EPSILON_SVR || model->param.svm_type == NU_SVR) && - model->probA!=NULL) - return model->probA[0]; - else - { - fprintf(stderr,"Model doesn't contain information for SVR probability inference\n"); - return 0; - } -} - -double svm_predict_values(const svm_model *model, const svm_node *x, double* dec_values) -{ - int i; - if(model->param.svm_type == ONE_CLASS || - model->param.svm_type == EPSILON_SVR || - model->param.svm_type == NU_SVR) - { - double *sv_coef = model->sv_coef[0]; - double sum = 0; - #pragma omp parallel for private(i) reduction(+:sum) schedule(guided) - for(i=0;il;i++) - sum += sv_coef[i] * Kernel::k_function(x,model->SV[i],model->param); - sum -= model->rho[0]; - *dec_values = sum; - - if(model->param.svm_type == ONE_CLASS) - return (sum>0)?1:-1; - else - return sum; - } - else - { - int nr_class = model->nr_class; - int l = model->l; - - double *kvalue = Malloc(double,l); - #pragma omp parallel for private(i) schedule(guided) - for(i=0;iSV[i],model->param); - - int *start = Malloc(int,nr_class); - start[0] = 0; - for(i=1;inSV[i-1]; - - int *vote = Malloc(int,nr_class); - for(i=0;inSV[i]; - int cj = model->nSV[j]; - - int k; - double *coef1 = model->sv_coef[j-1]; - double *coef2 = model->sv_coef[i]; - for(k=0;krho[p]; - dec_values[p] = sum; - - if(dec_values[p] > 0) - ++vote[i]; - else - ++vote[j]; - p++; - } - - int vote_max_idx = 0; - for(i=1;i vote[vote_max_idx]) - vote_max_idx = i; - - free(kvalue); - free(start); - free(vote); - return model->label[vote_max_idx]; - } -} - -double svm_predict(const svm_model *model, const svm_node *x) -{ - int nr_class = model->nr_class; - double *dec_values; - if(model->param.svm_type == ONE_CLASS || - model->param.svm_type == EPSILON_SVR || - model->param.svm_type == NU_SVR) - dec_values = Malloc(double, 1); - else - dec_values = Malloc(double, nr_class*(nr_class-1)/2); - double pred_result = svm_predict_values(model, x, dec_values); - free(dec_values); - return pred_result; -} - -double svm_predict_probability( - const svm_model *model, const svm_node *x, double *prob_estimates) -{ - if ((model->param.svm_type == C_SVC || model->param.svm_type == NU_SVC) && - model->probA!=NULL && model->probB!=NULL) - { - int i; - int nr_class = model->nr_class; - double *dec_values = Malloc(double, nr_class*(nr_class-1)/2); - svm_predict_values(model, x, dec_values); - - double min_prob=1e-7; - double **pairwise_prob=Malloc(double *,nr_class); - for(i=0;iprobA[k],model->probB[k]),min_prob),1-min_prob); - pairwise_prob[j][i]=1-pairwise_prob[i][j]; - k++; - } - if (nr_class == 2) - { - prob_estimates[0] = pairwise_prob[0][1]; - prob_estimates[1] = pairwise_prob[1][0]; - } - else - multiclass_probability(nr_class,pairwise_prob,prob_estimates); - - int prob_max_idx = 0; - for(i=1;i prob_estimates[prob_max_idx]) - prob_max_idx = i; - for(i=0;ilabel[prob_max_idx]; - } - else - return svm_predict(model, x); -} - -static const char *svm_type_table[] = -{ - "c_svc","nu_svc","one_class","epsilon_svr","nu_svr",NULL -}; - -static const char *kernel_type_table[]= -{ - "linear","polynomial","rbf","sigmoid","precomputed",NULL -}; - -int svm_save_model(const char *model_file_name, const svm_model *model) -{ - FILE *fp = fopen(model_file_name,"w"); - if(fp==NULL) return -1; - - char *old_locale = setlocale(LC_ALL, NULL); - if (old_locale) { - old_locale = strdup(old_locale); - } - setlocale(LC_ALL, "C"); - - const svm_parameter& param = model->param; - - fprintf(fp,"svm_type %s\n", svm_type_table[param.svm_type]); - fprintf(fp,"kernel_type %s\n", kernel_type_table[param.kernel_type]); - - if(param.kernel_type == POLY) - fprintf(fp,"degree %d\n", param.degree); - - if(param.kernel_type == POLY || param.kernel_type == RBF || param.kernel_type == SIGMOID) - fprintf(fp,"gamma %g\n", param.gamma); - - if(param.kernel_type == POLY || param.kernel_type == SIGMOID) - fprintf(fp,"coef0 %g\n", param.coef0); - - int nr_class = model->nr_class; - int l = model->l; - fprintf(fp, "nr_class %d\n", nr_class); - fprintf(fp, "total_sv %d\n",l); - - { - fprintf(fp, "rho"); - for(int i=0;irho[i]); - fprintf(fp, "\n"); - } - - if(model->label) - { - fprintf(fp, "label"); - for(int i=0;ilabel[i]); - fprintf(fp, "\n"); - } - - if(model->probA) // regression has probA only - { - fprintf(fp, "probA"); - for(int i=0;iprobA[i]); - fprintf(fp, "\n"); - } - if(model->probB) - { - fprintf(fp, "probB"); - for(int i=0;iprobB[i]); - fprintf(fp, "\n"); - } - - if(model->nSV) - { - fprintf(fp, "nr_sv"); - for(int i=0;inSV[i]); - fprintf(fp, "\n"); - } - - fprintf(fp, "SV\n"); - const double * const *sv_coef = model->sv_coef; - const svm_node * const *SV = model->SV; - - for(int i=0;ivalue)); - else - while(p->index != -1) - { - fprintf(fp,"%d:%.8g ",p->index,p->value); - p++; - } - fprintf(fp, "\n"); - } - - setlocale(LC_ALL, old_locale); - free(old_locale); - - if (ferror(fp) != 0 || fclose(fp) != 0) return -1; - else return 0; -} - -static char *line = NULL; -static int max_line_len; - -static char* readline(FILE *input) -{ - int len; - - if(fgets(line,max_line_len,input) == NULL) - return NULL; - - while(strrchr(line,'\n') == NULL) - { - max_line_len *= 2; - line = (char *) realloc(line,max_line_len); - len = (int) strlen(line); - if(fgets(line+len,max_line_len-len,input) == NULL) - break; - } - return line; -} - -// -// FSCANF helps to handle fscanf failures. -// Its do-while block avoids the ambiguity when -// if (...) -// FSCANF(); -// is used -// -#define FSCANF(_stream, _format, _var) do{ if (fscanf(_stream, _format, _var) != 1) return false; }while(0) -bool read_model_header(FILE *fp, svm_model* model) -{ - svm_parameter& param = model->param; - // parameters for training only won't be assigned, but arrays are assigned as NULL for safety - param.nr_weight = 0; - param.weight_label = NULL; - param.weight = NULL; - - char cmd[81]; - while(1) - { - FSCANF(fp,"%80s",cmd); - - if(strcmp(cmd,"svm_type")==0) - { - FSCANF(fp,"%80s",cmd); - int i; - for(i=0;svm_type_table[i];i++) - { - if(strcmp(svm_type_table[i],cmd)==0) - { - param.svm_type=i; - break; - } - } - if(svm_type_table[i] == NULL) - { - fprintf(stderr,"unknown svm type.\n"); - return false; - } - } - else if(strcmp(cmd,"kernel_type")==0) - { - FSCANF(fp,"%80s",cmd); - int i; - for(i=0;kernel_type_table[i];i++) - { - if(strcmp(kernel_type_table[i],cmd)==0) - { - param.kernel_type=i; - break; - } - } - if(kernel_type_table[i] == NULL) - { - fprintf(stderr,"unknown kernel function.\n"); - return false; - } - } - else if(strcmp(cmd,"degree")==0) - FSCANF(fp,"%d",¶m.degree); - else if(strcmp(cmd,"gamma")==0) - FSCANF(fp,"%lf",¶m.gamma); - else if(strcmp(cmd,"coef0")==0) - FSCANF(fp,"%lf",¶m.coef0); - else if(strcmp(cmd,"nr_class")==0) - FSCANF(fp,"%d",&model->nr_class); - else if(strcmp(cmd,"total_sv")==0) - FSCANF(fp,"%d",&model->l); - else if(strcmp(cmd,"rho")==0) - { - int n = model->nr_class * (model->nr_class-1)/2; - model->rho = Malloc(double,n); - for(int i=0;irho[i]); - } - else if(strcmp(cmd,"label")==0) - { - int n = model->nr_class; - model->label = Malloc(int,n); - for(int i=0;ilabel[i]); - } - else if(strcmp(cmd,"probA")==0) - { - int n = model->nr_class * (model->nr_class-1)/2; - model->probA = Malloc(double,n); - for(int i=0;iprobA[i]); - } - else if(strcmp(cmd,"probB")==0) - { - int n = model->nr_class * (model->nr_class-1)/2; - model->probB = Malloc(double,n); - for(int i=0;iprobB[i]); - } - else if(strcmp(cmd,"nr_sv")==0) - { - int n = model->nr_class; - model->nSV = Malloc(int,n); - for(int i=0;inSV[i]); - } - else if(strcmp(cmd,"SV")==0) - { - while(1) - { - int c = getc(fp); - if(c==EOF || c=='\n') break; - } - break; - } - else - { - fprintf(stderr,"unknown text in model file: [%s]\n",cmd); - return false; - } - } - - return true; - -} - -svm_model *svm_load_model(const char *model_file_name) -{ - FILE *fp = fopen(model_file_name,"rb"); - if(fp==NULL) return NULL; - - char *old_locale = setlocale(LC_ALL, NULL); - if (old_locale) { - old_locale = strdup(old_locale); - } - setlocale(LC_ALL, "C"); - - // read parameters - - svm_model *model = Malloc(svm_model,1); - model->rho = NULL; - model->probA = NULL; - model->probB = NULL; - model->sv_indices = NULL; - model->label = NULL; - model->nSV = NULL; - - // read header - if (!read_model_header(fp, model)) - { - fprintf(stderr, "ERROR: fscanf failed to read model\n"); - setlocale(LC_ALL, old_locale); - free(old_locale); - free(model->rho); - free(model->label); - free(model->nSV); - free(model); - return NULL; - } - - // read sv_coef and SV - - int elements = 0; - long pos = ftell(fp); - - max_line_len = 1024; - line = Malloc(char,max_line_len); - char *p,*endptr,*idx,*val; - - while(readline(fp)!=NULL) - { - p = strtok(line,":"); - while(1) - { - p = strtok(NULL,":"); - if(p == NULL) - break; - ++elements; - } - } - elements += model->l; - - fseek(fp,pos,SEEK_SET); - - int m = model->nr_class - 1; - int l = model->l; - model->sv_coef = Malloc(double *,m); - int i; - for(i=0;isv_coef[i] = Malloc(double,l); - model->SV = Malloc(svm_node*,l); - svm_node *x_space = NULL; - if(l>0) x_space = Malloc(svm_node,elements); - - int j=0; - for(i=0;iSV[i] = &x_space[j]; - - p = strtok(line, " \t"); - model->sv_coef[0][i] = strtod(p,&endptr); - for(int k=1;ksv_coef[k][i] = strtod(p,&endptr); - } - - while(1) - { - idx = strtok(NULL, ":"); - val = strtok(NULL, " \t"); - - if(val == NULL) - break; - x_space[j].index = (int) strtol(idx,&endptr,10); - x_space[j].value = strtod(val,&endptr); - - ++j; - } - x_space[j++].index = -1; - } - free(line); - - setlocale(LC_ALL, old_locale); - free(old_locale); - - if (ferror(fp) != 0 || fclose(fp) != 0) - return NULL; - - model->free_sv = 1; // XXX - return model; -} - -void svm_free_model_content(svm_model* model_ptr) -{ - if(model_ptr->free_sv && model_ptr->l > 0 && model_ptr->SV != NULL) - free((void *)(model_ptr->SV[0])); - if(model_ptr->sv_coef) - { - for(int i=0;inr_class-1;i++) - free(model_ptr->sv_coef[i]); - } - - free(model_ptr->SV); - model_ptr->SV = NULL; - - free(model_ptr->sv_coef); - model_ptr->sv_coef = NULL; - - free(model_ptr->rho); - model_ptr->rho = NULL; - - free(model_ptr->label); - model_ptr->label= NULL; - - free(model_ptr->probA); - model_ptr->probA = NULL; - - free(model_ptr->probB); - model_ptr->probB= NULL; - - free(model_ptr->sv_indices); - model_ptr->sv_indices = NULL; - - free(model_ptr->nSV); - model_ptr->nSV = NULL; -} - -void svm_free_and_destroy_model(svm_model** model_ptr_ptr) -{ - if(model_ptr_ptr != NULL && *model_ptr_ptr != NULL) - { - svm_free_model_content(*model_ptr_ptr); - free(*model_ptr_ptr); - *model_ptr_ptr = NULL; - } -} - -void svm_destroy_param(svm_parameter* param) -{ - free(param->weight_label); - free(param->weight); -} - -const char *svm_check_parameter(const svm_problem *prob, const svm_parameter *param) -{ - // svm_type - - int svm_type = param->svm_type; - if(svm_type != C_SVC && - svm_type != NU_SVC && - svm_type != ONE_CLASS && - svm_type != EPSILON_SVR && - svm_type != NU_SVR) - return "unknown svm type"; - - // kernel_type, degree - - int kernel_type = param->kernel_type; - if(kernel_type != LINEAR && - kernel_type != POLY && - kernel_type != RBF && - kernel_type != SIGMOID && - kernel_type != PRECOMPUTED) - return "unknown kernel type"; - - if(param->gamma < 0) - return "gamma < 0"; - - if(param->degree < 0) - return "degree of polynomial kernel < 0"; - - // cache_size,eps,C,nu,p,shrinking - - if(param->cache_size <= 0) - return "cache_size <= 0"; - - if(param->eps <= 0) - return "eps <= 0"; - - if(svm_type == C_SVC || - svm_type == EPSILON_SVR || - svm_type == NU_SVR) - if(param->C <= 0) - return "C <= 0"; - - if(svm_type == NU_SVC || - svm_type == ONE_CLASS || - svm_type == NU_SVR) - if(param->nu <= 0 || param->nu > 1) - return "nu <= 0 or nu > 1"; - - if(svm_type == EPSILON_SVR) - if(param->p < 0) - return "p < 0"; - - if(param->shrinking != 0 && - param->shrinking != 1) - return "shrinking != 0 and shrinking != 1"; - - if(param->probability != 0 && - param->probability != 1) - return "probability != 0 and probability != 1"; - - if(param->probability == 1 && - svm_type == ONE_CLASS) - return "one-class SVM probability output not supported yet"; - - - // check whether nu-svc is feasible - - if(svm_type == NU_SVC) - { - int l = prob->l; - int max_nr_class = 16; - int nr_class = 0; - int *label = Malloc(int,max_nr_class); - int *count = Malloc(int,max_nr_class); - - int i; - for(i=0;iy[i]; - int j; - for(j=0;jnu*(n1+n2)/2 > min(n1,n2)) - { - free(label); - free(count); - return "specified nu is infeasible"; - } - } - } - free(label); - free(count); - } - - return NULL; -} - -int svm_check_probability_model(const svm_model *model) -{ - return ((model->param.svm_type == C_SVC || model->param.svm_type == NU_SVC) && - model->probA!=NULL && model->probB!=NULL) || - ((model->param.svm_type == EPSILON_SVR || model->param.svm_type == NU_SVR) && - model->probA!=NULL); -} - -void svm_set_print_string_function(void (*print_func)(const char *)) -{ - if(print_func == NULL) - svm_print_string = &print_string_stdout; - else - svm_print_string = print_func; -} - -void svm_set_num_threads(int num_threads) -{ - omp_set_num_threads(num_threads); -} - -int svm_get_max_threads() -{ - return omp_get_max_threads(); -} diff --git a/deps/libsvm-3.22/svm.def b/deps/libsvm-3.22/svm.def deleted file mode 100644 index 0f25ecf..0000000 --- a/deps/libsvm-3.22/svm.def +++ /dev/null @@ -1,23 +0,0 @@ -LIBRARY libsvm -EXPORTS - svm_train @1 - svm_cross_validation @2 - svm_save_model @3 - svm_load_model @4 - svm_get_svm_type @5 - svm_get_nr_class @6 - svm_get_labels @7 - svm_get_svr_probability @8 - svm_predict_values @9 - svm_predict @10 - svm_predict_probability @11 - svm_free_model_content @12 - svm_free_and_destroy_model @13 - svm_destroy_param @14 - svm_check_parameter @15 - svm_check_probability_model @16 - svm_set_print_string_function @17 - svm_get_sv_indices @18 - svm_get_nr_sv @19 - svm_set_num_threads @20 - svm_get_max_threads @21 \ No newline at end of file diff --git a/deps/libsvm-3.22/svm.h b/deps/libsvm-3.22/svm.h deleted file mode 100644 index e05acea..0000000 --- a/deps/libsvm-3.22/svm.h +++ /dev/null @@ -1,108 +0,0 @@ -#ifndef _LIBSVM_H -#define _LIBSVM_H - -#define LIBSVM_VERSION 322 - -#ifdef __cplusplus -extern "C" { -#endif - -extern int libsvm_version; - -struct svm_node -{ - int index; - double value; -}; - -struct svm_problem -{ - int l; - double *y; - struct svm_node **x; -}; - -enum { C_SVC, NU_SVC, ONE_CLASS, EPSILON_SVR, NU_SVR }; /* svm_type */ -enum { LINEAR, POLY, RBF, SIGMOID, PRECOMPUTED }; /* kernel_type */ - -struct svm_parameter -{ - int svm_type; - int kernel_type; - int degree; /* for poly */ - double gamma; /* for poly/rbf/sigmoid */ - double coef0; /* for poly/sigmoid */ - - /* these are for training only */ - double cache_size; /* in MB */ - double eps; /* stopping criteria */ - double C; /* for C_SVC, EPSILON_SVR and NU_SVR */ - int nr_weight; /* for C_SVC */ - int *weight_label; /* for C_SVC */ - double* weight; /* for C_SVC */ - double nu; /* for NU_SVC, ONE_CLASS, and NU_SVR */ - double p; /* for EPSILON_SVR */ - int shrinking; /* use the shrinking heuristics */ - int probability; /* do probability estimates */ -}; - -// -// svm_model -// -struct svm_model -{ - struct svm_parameter param; /* parameter */ - int nr_class; /* number of classes, = 2 in regression/one class svm */ - int l; /* total #SV */ - struct svm_node **SV; /* SVs (SV[l]) */ - double **sv_coef; /* coefficients for SVs in decision functions (sv_coef[k-1][l]) */ - double *rho; /* constants in decision functions (rho[k*(k-1)/2]) */ - double *probA; /* pariwise probability information */ - double *probB; - int *sv_indices; /* sv_indices[0,...,nSV-1] are values in [1,...,num_traning_data] to indicate SVs in the training set */ - - /* for classification only */ - - int *label; /* label of each class (label[k]) */ - int *nSV; /* number of SVs for each class (nSV[k]) */ - /* nSV[0] + nSV[1] + ... + nSV[k-1] = l */ - /* XXX */ - int free_sv; /* 1 if svm_model is created by svm_load_model*/ - /* 0 if svm_model is created by svm_train */ -}; - -struct svm_model *svm_train(const struct svm_problem *prob, const struct svm_parameter *param); -void svm_cross_validation(const struct svm_problem *prob, const struct svm_parameter *param, int nr_fold, double *target); - -int svm_save_model(const char *model_file_name, const struct svm_model *model); -struct svm_model *svm_load_model(const char *model_file_name); - -int svm_get_svm_type(const struct svm_model *model); -int svm_get_nr_class(const struct svm_model *model); -void svm_get_labels(const struct svm_model *model, int *label); -void svm_get_sv_indices(const struct svm_model *model, int *sv_indices); -int svm_get_nr_sv(const struct svm_model *model); -double svm_get_svr_probability(const struct svm_model *model); - -double svm_predict_values(const struct svm_model *model, const struct svm_node *x, double* dec_values); -double svm_predict(const struct svm_model *model, const struct svm_node *x); -double svm_predict_probability(const struct svm_model *model, const struct svm_node *x, double* prob_estimates); - -void svm_free_model_content(struct svm_model *model_ptr); -void svm_free_and_destroy_model(struct svm_model **model_ptr_ptr); -void svm_destroy_param(struct svm_parameter *param); - -const char *svm_check_parameter(const struct svm_problem *prob, const struct svm_parameter *param); -int svm_check_probability_model(const struct svm_model *model); - -void svm_set_print_string_function(void (*print_func)(const char *)); - -void svm_set_num_threads(int num_threads); - -int svm_get_max_threads(); - -#ifdef __cplusplus -} -#endif - -#endif /* _LIBSVM_H */ diff --git a/src/LIBSVM.jl b/src/LIBSVM.jl index 7f14732..83876ec 100644 --- a/src/LIBSVM.jl +++ b/src/LIBSVM.jl @@ -1,9 +1,8 @@ __precompile__() module LIBSVM import LIBLINEAR -using Compat using SparseArrays -using Libdl +using libsvm_jll export svmtrain, svmpredict, fit!, predict, transform, SVC, NuSVC, OneClassSVM, NuSVR, EpsilonSVR, LinearSVC, @@ -44,7 +43,7 @@ end struct SVM{T} SVMtype::Type kernel::Kernel.KERNEL - weights::Union{Dict{T, Float64}, Compat.Nothing} + weights::Union{Dict{T, Float64}, Cvoid} nfeatures::Int nclasses::Int32 labels::Vector{T} @@ -158,45 +157,12 @@ function svmprint(str::Ptr{UInt8}) nothing end - -let libsvm = C_NULL - global get_libsvm - function get_libsvm() - if libsvm == C_NULL - if Sys.iswindows() - libsvm = Libdl.dlopen(joinpath(dirname(@__FILE__), "../deps", - "libsvm.dll")) - else - libsvm = Libdl.dlopen(joinpath(dirname(@__FILE__), "../deps", - "libsvm.so.2")) - end - ccall(Libdl.dlsym(libsvm, :svm_set_print_string_function), Compat.Nothing, - (Ptr{Compat.Nothing},), @cfunction(svmprint, Compat.Nothing, (Ptr{UInt8},) )) - end - libsvm - end -end - -macro cachedsym(symname) - cached = gensym() - quote - let $cached = C_NULL - global ($symname) - ($symname)() = ($cached) == C_NULL ? - ($cached = Libdl.dlsym(get_libsvm(), $(string(symname)))) : $cached - end - end +function __init__() + ccall((:svm_set_print_string_function, libsvm), Cvoid, + (Ptr{Cvoid},), @cfunction(svmprint, Cvoid, (Ptr{UInt8},) )) end -@cachedsym svm_train -@cachedsym svm_predict -@cachedsym svm_predict_values -@cachedsym svm_predict_probability -@cachedsym svm_free_model_content -@cachedsym svm_set_num_threads -@cachedsym svm_get_max_threads - function grp2idx(::Type{S}, labels::AbstractVector, label_dict::Dict{T, Int32}, reverse_labels::Vector{T}) where {T, S <: Real} @@ -257,7 +223,7 @@ end function indices_and_weights(labels::AbstractVector{T}, instances::AbstractMatrix{U}, - weights::Union{Dict{T, Float64}, Compat.Nothing}=nothing) where {T, U<:Real} + weights::Union{Dict{T, Float64}, Cvoid}=nothing) where {T, U<:Real} label_dict = Dict{T, Int32}() reverse_labels = Array{T}(undef, 0) idx = grp2idx(Float64, labels, label_dict, reverse_labels) @@ -293,10 +259,10 @@ function set_num_threads(nt::Integer) end if nt < 0 - nt = ccall(svm_get_max_threads(), Cint, ()) + nt = ccall((:svm_get_max_threads, libsvm), Cint, ()) end - ccall(svm_set_num_threads(), Compat.Nothing, (Cint,), nt) + ccall((:svm_set_num_threads, libsvm), Cvoid, (Cint,), nt) end """ @@ -306,7 +272,7 @@ svmtrain{T, U<:Real}(X::AbstractMatrix{U}, y::AbstractVector{T}=[]; gamma::Float64=1.0/size(X, 1), coef0::Float64=0.0, cost::Float64=1.0, nu::Float64=0.5, epsilon::Float64=0.1, tolerance::Float64=0.001, shrinking::Bool=true, - probability::Bool=false, weights::Union{Dict{T, Float64}, Compat.Nothing}=nothing, + probability::Bool=false, weights::Union{Dict{T, Float64}, Cvoid}=nothing, cachesize::Float64=200.0, verbose::Bool=false) ``` Train Support Vector Machine using LIBSVM using response vector `y` @@ -329,7 +295,7 @@ For one-class SVM use only `X`. * `tolerance::Float64=0.001`: tolerance of termination criterion * `shrinking::Bool=true`: whether to use the shrinking heuristics * `probability::Bool=false`: whether to train a SVC or SVR model for probability estimates -* `weights::Union{Dict{T, Float64}, Compat.Nothing}=nothing`: dictionary of class weights +* `weights::Union{Dict{T, Float64}, Cvoid}=nothing`: dictionary of class weights * `cachesize::Float64=100.0`: cache memory size in MB * `verbose::Bool=false`: print training output from LIBSVM if true * `nt::Integer=0`: number of OpenMP cores to use, if 0 it is set to OMP_NUM_THREADS, if negative it is set to the max number of threads @@ -342,7 +308,7 @@ function svmtrain(X::AbstractMatrix{U}, y::AbstractVector{T} = []; degree::Integer=3, gamma::Float64=1.0/size(X, 1), coef0::Float64=0.0, cost::Float64=1.0, nu::Float64=0.5, epsilon::Float64=0.1, tolerance::Float64=0.001, shrinking::Bool=true, - probability::Bool=false, weights::Union{Dict{T, Float64}, Compat.Nothing}=nothing, + probability::Bool=false, weights::Union{Dict{T, Float64}, Cvoid}=nothing, cachesize::Float64=200.0, verbose::Bool=false, nt::Integer=1) where {T, U<:Real} global verbosity @@ -381,12 +347,12 @@ function svmtrain(X::AbstractMatrix{U}, y::AbstractVector{T} = []; pointer(nodeptrs))] verbosity = verbose - mod = ccall(svm_train(), Ptr{SVMModel}, (Ptr{SVMProblem}, + mod = ccall((:svm_train, libsvm), Ptr{SVMModel}, (Ptr{SVMProblem}, Ptr{SVMParameter}), problem, param) svm = SVM(unsafe_load(mod), y, X, wts, reverse_labels, svmtype, kernel) - ccall(svm_free_model_content(), Compat.Nothing, (Ptr{Compat.Nothing},), mod) + ccall((:svm_free_model_content, libsvm), Cvoid, (Ptr{Cvoid},), mod) return (svm) #return(mod, weights, weight_labels) end @@ -426,15 +392,18 @@ function svmpredict(model::SVM{T}, X::AbstractMatrix{U}; nt::Integer=0) where {T end verbosity = false - fn = model.probability ? svm_predict_probability() : svm_predict_values() cmod, data = svmmodel(model) ma = [cmod] for i = 1:ninstances - output = ccall(fn, Float64, (Ptr{Compat.Nothing}, Ptr{SVMNode}, Ptr{Float64}), + if model.probability + output = ccall((:svm_predict_probability, libsvm ), Float64, (Ptr{Cvoid}, Ptr{SVMNode}, Ptr{Float64}), ma, nodeptrs[i], pointer(decvalues, nlabels*(i-1)+1)) - + else + output = ccall((:svm_predict_values, libsvm ), Float64, (Ptr{Cvoid}, Ptr{SVMNode}, Ptr{Float64}), + ma, nodeptrs[i], pointer(decvalues, nlabels*(i-1)+1)) + end if model.SVMtype == EpsilonSVR || model.SVMtype == NuSVR pred[i] = output elseif model.SVMtype == OneClassSVM diff --git a/src/ScikitLearnAPI.jl b/src/ScikitLearnAPI.jl index 0c12e2f..6c07edf 100644 --- a/src/ScikitLearnAPI.jl +++ b/src/ScikitLearnAPI.jl @@ -59,7 +59,7 @@ EpsilonSVR(;kernel::Kernel.KERNEL = Kernel.RadialBasis, gamma::Union{Float64,Sym @declare_hyperparameters(EpsilonSVR, [:kernel, :gamma, :epsilon, :cost, :degree, :coef0, :tolerance]) LinearSVC(;solver = Linearsolver.L2R_L2LOSS_SVC_DUAL, - weights::Union{Dict, Compat.Nothing} = nothing, tolerance::Float64=Inf, + weights::Union{Dict, Cvoid} = nothing, tolerance::Float64=Inf, cost::Float64 = 1.0, p::Float64 = 0.1, bias::Float64 = -1.0, verbose::Bool = false) = LinearSVC(solver, weights, tolerance, cost, p, bias, verbose, nothing) diff --git a/src/ScikitLearnTypes.jl b/src/ScikitLearnTypes.jl index 9b2e27f..3c52a3d 100644 --- a/src/ScikitLearnTypes.jl +++ b/src/ScikitLearnTypes.jl @@ -1,4 +1,4 @@ -using Compat + import ScikitLearnBase: BaseClassifier, BaseRegressor abstract type AbstractSVC<:BaseClassifier end @@ -7,7 +7,7 @@ abstract type AbstractSVR<:BaseRegressor end mutable struct SVC<:AbstractSVC kernel::Kernel.KERNEL gamma::Union{Float64,Symbol} - weights::Union{Dict, Compat.Nothing} + weights::Union{Dict, Cvoid} cost::Float64 degree::Int32 coef0::Float64 @@ -16,13 +16,13 @@ mutable struct SVC<:AbstractSVC probability::Bool verbose::Bool - fit::Union{SVM, Compat.Nothing} + fit::Union{SVM, Cvoid} end mutable struct NuSVC<:AbstractSVC kernel::Kernel.KERNEL gamma::Union{Float64,Symbol} - weights::Union{Dict, Compat.Nothing} + weights::Union{Dict, Cvoid} nu::Float64 cost::Float64 degree::Int32 @@ -31,7 +31,7 @@ mutable struct NuSVC<:AbstractSVC shrinking::Bool verbose::Bool - fit::Union{SVM, Compat.Nothing} + fit::Union{SVM, Cvoid} end mutable struct OneClassSVM<:AbstractSVC @@ -45,7 +45,7 @@ mutable struct OneClassSVM<:AbstractSVC shrinking::Bool verbose::Bool - fit::Union{SVM, Compat.Nothing} + fit::Union{SVM, Cvoid} end mutable struct NuSVR<:AbstractSVR @@ -59,7 +59,7 @@ mutable struct NuSVR<:AbstractSVR shrinking::Bool verbose::Bool - fit::Union{SVM, Compat.Nothing} + fit::Union{SVM, Cvoid} end mutable struct EpsilonSVR<:AbstractSVR @@ -73,7 +73,7 @@ mutable struct EpsilonSVR<:AbstractSVR shrinking::Bool verbose::Bool - fit::Union{SVM, Compat.Nothing} + fit::Union{SVM, Cvoid} end """ @@ -81,14 +81,14 @@ Linear SVM using LIBLINEAR """ mutable struct LinearSVC<:BaseClassifier solver::Linearsolver.LINEARSOLVER - weights::Union{Dict, Compat.Nothing} + weights::Union{Dict, Cvoid} tolerance::Float64 cost::Float64 p::Float64 bias::Float64 verbose::Bool - fit::Union{LIBLINEAR.LinearModel, Compat.Nothing} + fit::Union{LIBLINEAR.LinearModel, Cvoid} end #Map types to Int for Libsvm C api