Changeset 67


Ignore:
Timestamp:
Apr 27, 2004, 5:04:31 PM (19 years ago)
Author:
Peter
Message:

added choice to use weights for the calculation of kernel matrix

Location:
trunk/src
Files:
7 edited

Legend:

Unmodified
Added
Removed

  • trunk/src/GaussianKernelFunction.cc

    r42 r67  
    1717}
    1818
    19 double GaussianKernelFunction::operator()(const gslapi::vector& a,
    20                                           const gslapi::vector& b) const   
     19double GaussianKernelFunction::operator()(const gslapi::vector& a1,
     20                                          const gslapi::vector& a2,
     21                                          const gslapi::vector& w1,
     22                                          const gslapi::vector& w2) const   
    2123{
    22   u_int nof_valid_datapairs = 0;
    23   double tmp = 0;
    24   // check that a size == b size !
    25   if(a.size() != b.size())
    26     std::cerr << "vectors must have same length\n";
    27   for(u_int k=0; k<a.size(); k++)
    28     // Check that both a and b are Numbers (NaN is NOT equal to NaN)
    29     if (a(k) == a(k) & b(k) == b(k)) {
    30      tmp += (a(k)-b(k))*(a(k)-b(k));
    31      nof_valid_datapairs++;
    32     }
    33   //cout << nof_valid_datapairs << endl;
    34   if (nof_valid_datapairs == 0)
    35     std::cerr << "data file error: two columns "
    36               << " have no pair of valid data" << std::endl;
    37  
    38   tmp/=nof_valid_datapairs / a.size(); //correcting for missing values
    39  
    40  
    41   return exp(tmp/pow(sigma_,2));
     24  double distance = 0;
     25  double normalization_factor = 0;
     26  for (size_t i=0; i<a1.size(); i++) {
     27    distance += w1(i) * w2(i) * (a1(i)-a2(i)) * (a1(i)-a2(i));
     28    normalization_factor += w1(i) * w2(i);
     29  }
     30// to make it coherent with no weight case
     31  normalization_factor /= a1.size();
     32  return exp(distance/normalization_factor/pow(sigma_,2));
    4233}
    4334

  • trunk/src/GaussianKernelFunction.h

    r42 r67  
    55
    66#include "KernelFunction.h"
     7#include "vector.h"
    78
    89namespace theplu {
    … …  
    2425   
    2526  public:
    26     /**
    27        Constructor taking the sigma_ , i.e. the width of the Gaussian,
    28        as input. Default is sigma_ = 1.
    29     */
     27    ///Constructor taking the sigma_ , i.e. the width of the Gaussian,as
     28    ///input. Default is sigma_ = 1.
     29    ///
    3030    GaussianKernelFunction(double = 1);
    31     /**
    32        Destructor
    33     */
     31   
     32    ///Destructor
     33    ///
    3434       virtual ~GaussianKernelFunction(void) {};
    35     /**
    36        returning the scalar product of two vectors in feature space using
    37        the Gaussian kernel.
    38        @return
    39        \f$ exp((x - y)^{2}/\sigma^2) \f$ \n       
    40        
    41        
    42      */
    43     double operator()(const gslapi::vector&, const gslapi::vector&) const;
    4435   
     36    ///returning the scalar product of two vectors in feature space using the
     37    ///Gaussian kernel. @return \f$ exp((x - y)^{2}/\sigma^2) \f$ \n
     38    ///
     39    inline double operator()(const gslapi::vector& a1,
     40                             const gslapi::vector& a2)
     41      {return exp(-(a1-a2)*(a1-a2)/pow(sigma_,2));}
     42
     43    ///returning the scalar product of two vectors in feature space using the
     44    ///Gaussian kernel with weights.
     45    ///
     46    double operator()(const gslapi::vector& a1,
     47                      const gslapi::vector& a2,
     48                      const gslapi::vector& w1,
     49                      const gslapi::vector& w2) const;
     50         
    4551  private:
    4652    double sigma_;

  • trunk/src/Kernel.cc

    r51 r67  
    1313  : k_(data.rows(),data.rows())
    1414{
    15   for(u_int i=0;i<data.rows();i++)
    16     for(u_int j=0;j<i+1;j++)
     15  for(u_int i=0;i<data.rows();i++)
     16    for(u_int j=0;j<i+1;j++)
    1717      k_(i,j)=kf(data[i],data[j]);
    1818 
    … …  
    2323}
    2424
     25Kernel::Kernel(const gslapi::matrix& data, const KernelFunction& kf,
     26               const gslapi::matrix& weight)
     27  : k_(data.rows(),data.rows())
     28{
     29  for(u_int i=0;i<data.rows();i++)
     30    for(u_int j=0;j<i+1;j++)
     31      k_(i,j)=kf(data[i],data[j],weight[i],weight[j]);
     32 
     33  // Copy lower triangle to upper triangle of Kernel matrix
     34  for(u_int i=0;i<data.rows()-1;i++)
     35    for(u_int j=i+1;j<data.rows();j++)
     36      k_(i,j)=k_(j,i);
     37}
     38
     39
    2540Kernel::~Kernel(void)
    2641{

  • trunk/src/Kernel.h

    r51 r67  
    2222  public:
    2323    ///
    24     ///   Constructor taking the data matrix and KernelFunction as input.
     24    ///   Constructor taking the data matrix and KernelFunction as
     25    ///   input. Note: can not handle NaNs. When dealing with missing values,
     26    ///   use constructor taking a weight matrix.
     27    Kernel(const gslapi::matrix&, const KernelFunction&);
     28   
    2529    ///
    26     Kernel(const gslapi::matrix&, const KernelFunction&);
     30    ///   Constructor taking the data matrix, the KernelFunction and a weight
     31    ///   matrix as input.
     32    Kernel(const gslapi::matrix&, const KernelFunction&, const gslapi::matrix&);
     33
    2734    ///
    2835    ///   Destructor
    29     ///
    3036    virtual ~Kernel(void);
    31     ///
     37
    3238    ///   @return Kernel matrix
    3339    ///

  • trunk/src/KernelFunction.h

    r42 r67  
    3333                              const gslapi::vector&) const = 0;
    3434   
     35    virtual double operator()(const gslapi::vector&,
     36                              const gslapi::vector&,
     37                              const gslapi::vector&,
     38                              const gslapi::vector&) const = 0;
     39   
    3540  }; // class KernelFunction
    3641

  • trunk/src/PolynomialKernelFunction.cc

    r42 r67  
    1717}
    1818
    19 double PolynomialKernelFunction::operator()(const gslapi::vector& a,
    20                                             const gslapi::vector& b) const   
     19double PolynomialKernelFunction::operator()(const gslapi::vector& a1,
     20                                            const gslapi::vector& a2) const   
    2121{
    22   using namespace std;
     22  if(order_>1)
     23    return pow(a1*a2,order_);
     24  return a1*a2;
     25}
    2326
    24   u_int nof_valid_datapairs = 0;
    25   double tmp = 0;
    26   // check that a size == b size !
    27   if(a.size() != b.size())
    28     cerr << "vectors must have same length\n";
    29   for(u_int k=0; k<a.size(); k++)
    30     // Check that both a and b are Numbers (NaN is NOT equal to NaN)
    31     if (a(k) == a(k) & b(k) == b(k)) {
    32      tmp += a(k)*b(k);
    33      nof_valid_datapairs++;
    34     }
    35   //cout << nof_valid_datapairs << endl;
    36   if (nof_valid_datapairs == 0)
    37     cerr << "data file error: two columns "
    38    << " have no pair of valid data\n";
    39  
    40   tmp/=nof_valid_datapairs / a.size(); //correcting for missing values
    41  
    42   if(order_>1)
    43     return pow(1+tmp,order_);
    44   return tmp;
     27double PolynomialKernelFunction::operator()(const gslapi::vector& a1,
     28                                            const gslapi::vector& a2,
     29                                            const gslapi::vector& w1,
     30                                            const gslapi::vector& w2) const   
     31{
     32  double scalar = 0;
     33  double normalization_factor = 0;
     34  for (size_t i=0; i<a1.size(); i++) {
     35    scalar += w1(i) * w2(i) * a1(i) * a2(i);
     36    normalization_factor += w1(i) * w2(i);
     37  }
     38// to make it coherent with no weight case
     39  normalization_factor /= a1.size();
     40  if(order_>1)
     41    return pow(scalar/normalization_factor,order_);
     42  return scalar/normalization_factor;
    4543}
    4644

  • trunk/src/PolynomialKernelFunction.h

    r42 r67  
    1414namespace cpptools {
    1515
    16   /**
    17      Class calculating one element in the kernel matrix using the
    18      polynomial kernel.
    19   */
     16  ///Class calculating one element in the kernel matrix using the polynomial
     17  ///kernel.
     18  ///
    2019
    2120 
    … …  
    2423   
    2524  public:
    26     /**
    27        Constructor taking the order of the polynomial as input. Default
    28        is order=1 yielding the linear kernel function.
    29     */
     25    ///Constructor taking the order of the polynomial as input. Default is
     26    ///order=1 yielding the linear kernel function.
     27    ///
    3028    PolynomialKernelFunction(int = 1);
    31     /**
    32        Destructor
    33     */
     29   
     30    ///Destructor
     31    ///
    3432       virtual ~PolynomialKernelFunction(void) {};
    35     /**
    36        returning the scalar product of two vectors in feature space using
    37        the polynomial kernel.
    38        @return
    39        If order is larger than one: \f$ (1+x \cdot y)^{order} \f$ \n       
    40        If order is one (linear):    \f$ x \cdot y \f$ 
    41        
    42      */
     33   
     34    ///returning the scalar product of two vectors in feature space using the
     35    ///polynomial kernel. @return If order is larger than one: \f$ (1+x \cdot
     36    ///y)^{order} \f$ \n If order is one (linear): \f$ x \cdot y \f$
     37    ///   
    4338    double operator()(const gslapi::vector&, const gslapi::vector&) const;
     39
     40    ///returning the scalar product of two vectors in feature space using the
     41    ///polynomial kernel with weights.
     42    ///   
     43    double operator()(const gslapi::vector&, const gslapi::vector&,
     44                      const gslapi::vector&, const gslapi::vector&) const;
    4445   
    4546  private:
Note: See TracChangeset for help on using the changeset viewer.