source: trunk/src/SVM.cc @ 30

// $Id: SVM.cc 30 2004-01-16 17:40:51Z peter $

#ifndef _THEP_SVM_
#define _THEP_SVM_

// System includes
#include <math.h>

// Thep C++ Tools
#include "SVM.h"
#include "matrix.h"
#include "vector.h"
#include "random_singleton.h"

using namespace thep_cpp_tools;
using namespace std;

SVM::SVM(const thep_gsl_api::matrix& kernel, 
   const thep_gsl_api::vector& target) : trained_(false), 
                 kernel_(kernel),
                 target_(target),
                 alpha_(target.size(),true,true)
{
}

void SVM::train() //Should be done so one can choose to not train on all the samples
{
  thep_gsl_api::vector E = thep_gsl_api::vector(-target_);
  
  double upper_bound = pow(10., 32); 
  u_int count = 0;
  double alpha_new;
  random_singleton* rnd=random_singleton::get_instance();
  double u;
  double v;
  
  // Stop criteria  
  ofstream myout("alpha.tst");
  double tmp;
  thep_gsl_api::vector dalpha;  
  thep_gsl_api::vector one(E.size(),true,true);
  one.set_all(1);
  while (count<10000)
    {
     count++;  
     dalpha = alpha_.mul_elements(target_);
     E=kernel_*dalpha-target_; //should be done in another way!!
     tmp=dalpha*(kernel_*dalpha)/2;
     tmp=tmp-alpha_*one;
     myout << tmp << '\n';
     // Choosing a pair of variables to modify (Should be done more clever)
     u_long index1 = rnd->get_uniform_int(kernel_.cols());
     u_long index2 = rnd->get_uniform_int(kernel_.cols()-1);
     if (index2 >= index1) 
      index2++;
  
     //Updating the two variables
     if (target_[index1]!=target_[index2])
       {       
  if (alpha_[index2] > alpha_[index1] )
    {
    v = upper_bound;
    u = alpha_[index2] - alpha_[index1];
    }
  else
    {
     v = upper_bound - alpha_[index1] + alpha_[index2];
     u = 0;
    }
       }
     else
       {       
  if (alpha_[index2] + alpha_[index1] > upper_bound)
    {
    u = alpha_[index2] + alpha_[index1] - upper_bound;
    v = upper_bound;    
    }
  else
    {
     u = 0;
     v = alpha_[index1] + alpha_[index2];
    }
       }

     double k = kernel_.get(index1, index1) + kernel_.get(index2, index2) - 
       2*kernel_.get(index1, index2);
     alpha_new = alpha_[index2] + target_[index2]*
       (E[index1]-E[index2])/k;
     if (alpha_new > v)
       alpha_new = v;
     else if (alpha_new<u)
       alpha_new = u;
     
     alpha_[index1]+=target_[index1]*target_[index2]*(alpha_[index2]-alpha_new);
     alpha_[index2]=alpha_new;
     
     
     
     
    

    }
  
  myout.close();
  
  //thep_gsl_api::vector output = kernel_*dalpha;
  //for (int i=0; i<64; i++){ 
  //cout << output.get(i) << '\t';
  //cout << alpha_.get(i) << endl;
  //}

  trained_= true;
}

#endif