added Matlab interface, removed IL options for plain GPHIK

FMKGPHyperparameterOptimization.h

@@ -1,6 +1,6 @@
 /** 
 * @file FMKGPHyperparameterOptimization.h
-* @brief Heart of the framework to set up everything, perform optimization, incremental updates, classification, variance prediction (Interface)
+* @brief Heart of the framework to set up everything, perform optimization, classification, and variance prediction (Interface)
 * @author Erik Rodner, Alexander Freytag
 * @date 01/02/2012
 
@@ -33,7 +33,7 @@ namespace NICE {
   
   /** 
  * @class FMKGPHyperparameterOptimization
- * @brief Heart of the framework to set up everything, perform optimization, incremental updates, classification, variance prediction
+ * @brief Heart of the framework to set up everything, perform optimization, classification, and variance prediction
  * @author Erik Rodner, Alexander Freytag
  */
   
@@ -102,52 +102,40 @@ class FMKGPHyperparameterOptimization : NICE::Persistent
     double * precomputedTForVarEst;
 
     //! optimize noise with the GP likelihood
-    bool optimizeNoise;
-
-    //! learn in a balanced manner
-    bool learnBalanced;       
+    bool optimizeNoise;     
        
-    //! k largest eigenvalues for every kernel matrix (k == nrOfEigenvaluesToConsider, if we do not use balanced learning, we have only 1 entry)
-    std::vector< NICE::Vector> eigenMax;
+    //! k largest eigenvalues of the kernel matrix (k == nrOfEigenvaluesToConsider)
+    NICE::Vector eigenMax;
 
-    //! eigenvectors corresponding to k largest eigenvalues for every matrix (k == nrOfEigenvaluesToConsider) -- format: nxk
-    std::vector< NICE::Matrix> eigenMaxVectors;
+    //! eigenvectors corresponding to k largest eigenvalues (k == nrOfEigenvaluesToConsider) -- format: nxk
+    NICE::Matrix eigenMaxVectors;
     
     //! needed for optimization and variance approximation
-    std::map<int, IKMLinearCombination * > ikmsums;
+    IKMLinearCombination * ikmsum;
     
     //! storing the labels is needed for Incremental Learning (re-optimization)
     NICE::Vector labels;
     
-    //! we store the alpha vectors for good initializations in the IL setting
-    std::map<int, NICE::Vector> lastAlphas;
 
     //! calculate binary label vectors using a multi-class label vector
     int prepareBinaryLabels ( std::map<int, NICE::Vector> & binaryLabels, const NICE::Vector & y , std::set<int> & myClasses);     
     
-    //! prepare the GPLike objects for given binary labels and already given ikmsum-objects
-    inline void setupGPLikelihoodApprox(std::map<int,GPLikelihoodApprox * > & gplikes, const std::map<int, NICE::Vector> & binaryLabels, std::map<int,uint> & parameterVectorSizes);    
+    //! prepare the GPLike object for given binary labels and already given ikmsum-object
+    inline void setupGPLikelihoodApprox( GPLikelihoodApprox * & gplike, const std::map<int, NICE::Vector> & binaryLabels, uint & parameterVectorSize);    
     
     //! update eigenvectors and eigenvalues for given ikmsum-objects and a method to compute eigenvalues
-    inline void updateEigenVectors();
+    inline void updateEigenDecomposition( const int & i_noEigenValues );
     
     //! core of the optimize-functions
-    inline void performOptimization(std::map<int,GPLikelihoodApprox * > & gplikes, const std::map<int,uint> & parameterVectorSizes, const bool & roughOptimization = true);
+    inline void performOptimization( GPLikelihoodApprox & gplike, const uint & parameterVectorSize);
     
     //! apply the optimized transformation values to the underlying features
-    inline void transformFeaturesWithOptimalParameters(const std::map<int,GPLikelihoodApprox * > & gplikes, const std::map<int,uint> & parameterVectorSizes);
+    inline void transformFeaturesWithOptimalParameters(const GPLikelihoodApprox & gplike, const uint & parameterVectorSize);
     
     //! build the resulting matrices A and B as well as lookup tables T for fast evaluations using the optimized parameter settings
-    inline void computeMatricesAndLUTs(const std::map<int,GPLikelihoodApprox * > & gplikes);
-    
-    //! Update matrices (A, B, LUTs) and optionally find optimal parameters after adding a new example.  
-    void updateAfterSingleIncrement (const NICE::SparseVector & x, const bool & performOptimizationAfterIncrement = false);    
-    //! Update matrices (A, B, LUTs) and optionally find optimal parameters after adding multiple examples.  
-    void updateAfterMultipleIncrements (const std::vector<const NICE::SparseVector*> & x, const bool & performOptimizationAfterIncrement = false);   
-    
-    //! use the alphas from the last iteration as initial guess for the ILS?
-    bool usePreviousAlphas;
+    inline void computeMatricesAndLUTs( const GPLikelihoodApprox & gplike);
     
+     
     //! store the class number of the positive class (i.e., larger class no), only used in binary settings
     int binaryLabelPositive;
     //! store the class number of the negative class (i.e., smaller class no), only used in binary settings
@@ -155,8 +143,7 @@ class FMKGPHyperparameterOptimization : NICE::Persistent
     
     //! contains all class numbers of the currently known classes
     std::set<int> knownClasses;
-    //! contains the class numbers of new classes - only needed within the increment step
-    std::set<int> newClasses;
+
     
   public:  
     
@@ -179,6 +166,8 @@ class FMKGPHyperparameterOptimization : NICE::Persistent
     void setParameterUpperBound(const double & _parameterUpperBound);
     void setParameterLowerBound(const double & _parameterLowerBound);  
     
+    std::set<int> getKnownClassNumbers ( ) const;
+    
     //high level methods
     
     void initialize( const Config *conf, ParameterizedFunction *pf, FastMinKernel *fmk = NULL, const std::string & confSection = "GPHIKClassifier" );
@@ -219,11 +208,18 @@ class FMKGPHyperparameterOptimization : NICE::Persistent
     void optimize ( std::map<int, NICE::Vector> & binaryLabels );    
     
     /**
-    * @brief Compute the necessary variables for appxorimations of predictive variance, assuming an already initialized fmk object
+    * @brief Compute the necessary variables for appxorimations of predictive variance (LUTs), assuming an already initialized fmk object
     * @author Alexander Freytag
     * @date 11-04-2012 (dd-mm-yyyy)
     */       
-    void prepareVarianceApproximation();
+    void prepareVarianceApproximationRough();
+    
+    /**
+    * @brief Compute the necessary variables for fine appxorimations of predictive variance (EVs), assuming an already initialized fmk object
+    * @author Alexander Freytag
+    * @date 11-04-2012 (dd-mm-yyyy)
+    */       
+    void prepareVarianceApproximationFine();    
     
     /**
     * @brief classify an example 
@@ -249,43 +245,88 @@ class FMKGPHyperparameterOptimization : NICE::Persistent
     */
     int classify ( const NICE::Vector & x, SparseVector & scores ) const;    
 
+    //////////////////////////////////////////
+    // variance computation: sparse inputs
+    //////////////////////////////////////////
+    
     /**
     * @brief compute predictive variance for a given test example using a rough approximation: k_{**} -  k_*^T (K+\sigma I)^{-1} k_* <= k_{**} - |k_*|^2 * 1 / \lambda_max(K + \sigma I), where we approximate |k_*|^2 by neglecting the mixed terms
     * @author Alexander Freytag
     * @date 10-04-2012 (dd-mm-yyyy)
     * @param x input example
-    * @param predVariances contains the approximations of the predictive variances
+    * @param predVariance contains the approximation of the predictive variance
     *
     */    
-    void computePredictiveVarianceApproximateRough(const NICE::SparseVector & x, NICE::Vector & predVariances) const;
-
+    void computePredictiveVarianceApproximateRough(const NICE::SparseVector & x, double & predVariance ) const;
+    
     /**
     * @brief compute predictive variance for a given test example using a fine approximation  (k eigenvalues and eigenvectors to approximate the quadratic term)
     * @author Alexander Freytag
     * @date 18-04-2012 (dd-mm-yyyy)
     * @param x input example
-    * @param predVariances contains the approximations of the predictive variances
+     * @param predVariance contains the approximation of the predictive variance
     *
     */    
-    void computePredictiveVarianceApproximateFine(const NICE::SparseVector & x, NICE::Vector & predVariances) const;    
+    void computePredictiveVarianceApproximateFine(const NICE::SparseVector & x, double & predVariance ) const; 
     
     /**
     * @brief compute exact predictive variance for a given test example using ILS methods (exact, but more time consuming than approx versions)
     * @author Alexander Freytag
     * @date 10-04-2012 (dd-mm-yyyy)
     * @param x input example
-    * @param predVariances contains the approximations of the predictive variances
+     * @param predVariance contains the approximation of the predictive variance
+    *
+    */    
+    void computePredictiveVarianceExact(const NICE::SparseVector & x, double & predVariance ) const; 
+    
+    
+    //////////////////////////////////////////
+    // variance computation: non-sparse inputs
+    //////////////////////////////////////////
+    
+    /**
+    * @brief compute predictive variance for a given test example using a rough approximation: k_{**} -  k_*^T (K+\sigma I)^{-1} k_* <= k_{**} - |k_*|^2 * 1 / \lambda_max(K + \sigma I), where we approximate |k_*|^2 by neglecting the mixed terms
+    * @author Alexander Freytag
+    * @date 19-12-2013 (dd-mm-yyyy)
+    * @param x input example
+    * @param predVariance contains the approximation of the predictive variance
     *
     */    
-    void computePredictiveVarianceExact(const NICE::SparseVector & x, NICE::Vector & predVariances) const;
+    void computePredictiveVarianceApproximateRough(const NICE::Vector & x, double & predVariance ) const;    
+
+   
+    
+    /**
+    * @brief compute predictive variance for a given test example using a fine approximation  (k eigenvalues and eigenvectors to approximate the quadratic term)
+    * @author Alexander Freytag
+    * @date 19-12-2013 (dd-mm-yyyy)
+    * @param x input example
+     * @param predVariance contains the approximation of the predictive variance
+    *
+    */    
+    void computePredictiveVarianceApproximateFine(const NICE::Vector & x, double & predVariance ) const;      
+    
+
+    
+   /**
+    * @brief compute exact predictive variance for a given test example using ILS methods (exact, but more time consuming than approx versions)
+    * @author Alexander Freytag
+    * @date 19-12-2013 (dd-mm-yyyy)
+    * @param x input example
+    * @param predVariance contains the approximation of the predictive variance
+    *
+    */    
+    void computePredictiveVarianceExact(const NICE::Vector & x, double & predVariance ) const;  
+    
+    
+    
+    
     
     /** Persistent interface */
     void restore ( std::istream & is, int format = 0 );
     void store ( std::ostream & os, int format = 0 ) const;
     void clear ( ) ;
     
-    void addExample( const NICE::SparseVector & x, const double & label, const bool & performOptimizationAfterIncrement = true);
-    void addMultipleExamples( const std::vector<const NICE::SparseVector*> & newExamples, const NICE::Vector & labels, const bool & performOptimizationAfterIncrement = false);
         
 };
 

FastMinKernel.cpp

@@ -1142,6 +1142,10 @@ double* FastMinKernel::hikPrepareLookupTableForKVNApproximation(const Quantizati
   return Tlookup;  
 }
 
+    //////////////////////////////////////////
+    // variance computation: sparse inputs
+    //////////////////////////////////////////    
+
 void FastMinKernel::hikComputeKVNApproximation(const NICE::VVector & A, const NICE::SparseVector & xstar, double & norm, const ParameterizedFunction *pf ) 
 {
   norm = 0.0;
@@ -1227,6 +1231,121 @@ void FastMinKernel::hikComputeKernelVector ( const NICE::SparseVector& xstar, NI
     }
     
 
+    int position;
+
+    //where is the example x^z_i located in
+    //the sorted array? -> perform binary search, runtime O(log(n))
+    // search using the original value
+    X_sorted.findFirstLargerInDimension(dim, fval, position);
+    position--;
+    
+    //get the non-zero elements for this dimension  
+    const multimap< double, SortedVectorSparse<double>::dataelement> & nonzeroElements = X_sorted.getFeatureValues(dim).nonzeroElements();
+    
+    //run over the non-zero elements and add the corresponding entries to our kernel vector
+
+    int count(nrZeroIndices);
+    for ( SortedVectorSparse<double>::const_elementpointer i = nonzeroElements.begin(); i != nonzeroElements.end(); i++, count++ )
+    {
+      int origIndex(i->second.first); //orig index (i->second.second would be the transformed feature value)
+      if (count <= position)
+        kstar[origIndex] += i->first; //orig feature value
+      else
+        kstar[origIndex] += fval;
+    }
+  }  
+}
+
+    //////////////////////////////////////////
+    // variance computation: non-sparse inputs
+    //////////////////////////////////////////  
+
+void FastMinKernel::hikComputeKVNApproximation(const NICE::VVector & A, const NICE::Vector & xstar, double & norm, const ParameterizedFunction *pf ) 
+{
+  norm = 0.0;
+  int dim ( 0 );
+  for (Vector::const_iterator i = xstar.begin(); i != xstar.end(); i++, dim++)
+  {
+  
+    double fval = *i;
+    
+    int nrZeroIndices = X_sorted.getNumberOfZeroElementsPerDimension(dim);
+    if ( nrZeroIndices == n ) {
+      // all features are zero so let us ignore them completely
+      continue;
+    }
+
+    int position;
+
+    //where is the example x^z_i located in
+    //the sorted array? -> perform binary search, runtime O(log(n))
+    // search using the original value
+    X_sorted.findFirstLargerInDimension(dim, fval, position);
+    position--;
+  
+    //NOTE again - pay attention! This is only valid if all entries are NOT negative! - if not, ask wether the current feature is greater than zero. If so, subtract the nrZeroIndices, if not do not
+    double firstPart(0.0);
+    //TODO in the "overnext" line there occurs the following error
+    // Invalid read of size 8    
+    if (position >= 0) 
+      firstPart = (A[dim][position-nrZeroIndices]);
+    else
+      firstPart = 0.0;
+    
+    double secondPart( 0.0);
+      
+    if ( pf != NULL )
+      fval = pf->f ( dim, fval );
+    
+    fval = fval * fval;
+    
+    if (position >= 0) 
+      secondPart = fval * (n-nrZeroIndices-(position+1));
+    else //if x_d^* is smaller than every non-zero training example
+      secondPart = fval * (n-nrZeroIndices);
+    
+    // but apply using the transformed one
+    norm += firstPart + secondPart;
+  }  
+}
+
+void FastMinKernel::hikComputeKVNApproximationFast(const double *Tlookup, const Quantization & q, const NICE::Vector & xstar, double & norm) const
+{
+  norm = 0.0;
+  // runtime is O(d) if the quantizer is O(1)
+  int dim ( 0 );
+  for (Vector::const_iterator i = xstar.begin(); i != xstar.end(); i++, dim++ )
+  {
+    double v = *i;
+    // we do not need a parameterized function here, since the quantizer works on the original feature values. 
+    // nonetheless, the lookup table was created using the parameterized function    
+    uint qBin = q.quantize(v);
+    
+    norm += Tlookup[dim*q.size() + qBin];
+  }  
+}
+
+
+void FastMinKernel::hikComputeKernelVector( const NICE::Vector & xstar, NICE::Vector & kstar) const
+{
+  //init
+  kstar.resize(this->n);
+  kstar.set(0.0);
+  
+  //let's start :)
+  int dim ( 0 );
+  for (NICE::Vector::const_iterator i = xstar.begin(); i != xstar.end(); i++, dim++)
+  {
+  
+    double fval = *i;
+    
+    int nrZeroIndices = X_sorted.getNumberOfZeroElementsPerDimension(dim);
+    if ( nrZeroIndices == n ) {
+      // all features are zero so let us ignore them completely
+      continue;
+    }
+    
+
     int position;
 
     //where is the example x^z_i located in
@@ -1327,214 +1446,3 @@ bool FastMinKernel::getDebug( )   const
 {
   return debug;
 }
-
-// ----------------- INCREMENTAL LEARNING METHODS -----------------------
-void FastMinKernel::addExample(const NICE::SparseVector & _v, const ParameterizedFunction *pf )
-{
-  X_sorted.add_feature(_v, pf );
-  n++;
-}
-void FastMinKernel::addExample(const std::vector<double> & _v, const ParameterizedFunction *pf )
-{
-  X_sorted.add_feature(_v, pf );
-  n++;
-}
-
-void FastMinKernel::updatePreparationForAlphaMultiplications(const NICE::SparseVector & _v, const double & alpha, NICE::VVector & A, NICE::VVector & B, const ParameterizedFunction *pf) const
-{ 
-  NICE::SparseVector::const_iterator it = _v.begin();
-  for (int dim = 0; dim < this->d; dim++)
-  {
-    if (it->first == dim)
-    {
-      //increase both datastructures by one
-      A[dim].append(0.0);
-      B[dim].append(0.0);
-      
-      //this is the index of the new example in this dimension, which was already added
-      int idx;
-      X_sorted.findLastInDimension(dim, it->second, idx);
-      //actually we do not want to have the next position, but the current one
-      idx--;
-      
-      // and we do not care about zero elements since we store matrices A and B only for non-zero elements in the training data
-      idx -= X_sorted.getNumberOfZeroElementsPerDimension(dim);
-      
-      // we start using the last old element, which is located at size-2
-      for( int i = A[dim].size()-2; i >= (idx-1); i--)
-      {
-        if (pf != NULL)
-          A[dim][i+1] = A[dim][i] + alpha * pf->f ( 1, it->second );
-        else
-        {
-          A[dim][i+1] = A[dim][i] + alpha * it->second;
-        }
-      }    
-          
-      // remember: in contrast to the explanations in our ECCV-paper, we store the alpha-values of the INCREASINGLY ordered features
-      // in the matrix B, not in decreasing order
-      for (int i = B[dim].size()-1; i >= std::max(1,idx); i--)
-      {
-        B[dim][i] = B[dim][i-1] + alpha;
-      }
-      
-      //special case
-      if (idx == 0)
-      {
-        if (pf != NULL)
-          A[dim][0] = alpha * pf->f ( 1, it->second );
-        else
-          A[dim][0] = alpha * it->second;
-        
-        B[dim][0] = alpha;
-      }      
-      
-      it++;
-    }
-    else //_v is zero for that dimension
-    {
-      //nothing to do, since we do not store any information about zero elements
-    }
-  }
-}
-
-void FastMinKernel::updateLookupTableForAlphaMultiplications(const NICE::SparseVector & _v, const double & alpha, double * T, const Quantization & q, const ParameterizedFunction *pf) const
-{
-  //be aware, index n-1 is only valid, if we do not explicitely changed the indices while inserting elements
-  //actually, the code written below works equally to the following line, but is more efficient since we do not have to call the feature matrix several times
-//   this->hikUpdateLookupTable(T, alpha, 0.0, n-1, q, pf );
-  if (T == NULL)
-  {
-    fthrow(Exception, "FastMinKernel::updateLookupTableForAlphaMultiplications LUT not initialized, run FastMinKernel::hikPrepareLookupTable first!");
-    return;
-  }
-  
-  // number of quantization bins
-  uint hmax = q.size();
-
-  // store (transformed) prototypes
-  double *prototypes = new double [ hmax ];
-  for ( uint i = 0 ; i < hmax ; i++ )
-    if ( pf != NULL ) {
-      // FIXME: the transformed prototypes could change from dimension to another dimension
-      // We skip this flexibility ...but it should be changed in the future
-      prototypes[i] = pf->f ( 1, q.getPrototype(i) );
-    } else {
-      prototypes[i] = q.getPrototype(i);
-    }
-  
-  // loop through all dimensions
-  for (NICE::SparseVector::const_iterator it = _v.begin(); it != _v.end(); it++)
-  {
-
-    int dim(it->first);
-
-    double x_i = it->second;
-    //as usually, we quantize the original features, but use the quantized transformed features lateron
-    int q_bin = q.quantize(x_i);      
-
-    //TODO we could speed up this with first do a binary search for the position where the min changes, and then do two separate for-loops
-    for (uint j = 0; j < hmax; j++)
-    {
-      double fval;
-      
-      if (q_bin > j)
-        fval = prototypes[j]; //the prototypes are already transformed
-      else
-      {
-        if (pf != NULL)
-          fval = pf->f( 1, x_i );
-        else
-          fval = x_i;
-      }
-      
-      // pay attention: we use either the quantized prototypes or the REAL feature values, not the quantized ones!
-      T[ dim*hmax + j ] += alpha*fval;
-    }
-  }
-
-  delete [] prototypes;
-}
-
-void FastMinKernel::updatePreparationForKVNApproximation(const NICE::SparseVector & _v, NICE::VVector & A, const ParameterizedFunction *pf) const
-{
-  for (NICE::SparseVector::const_iterator it = _v.begin(); it != _v.end(); it++)
-  {
-    int dim(it->first);  
-    int idx;
-    
-    // we use the original feature value for this search, not the transformed one (see FeatureMatrixT)
-    // we assume that the nex example was already inserted to the FeatureMatrix
-    X_sorted.findLastInDimension(dim, it->second, idx);  
-    //we do not want to considere zero elements, since we store it in a sparse way
-    idx -= X_sorted.getNumberOfZeroElementsPerDimension(dim);
-    // not the next one, but the current (position vs index)
-    idx--;
-    
-    // perform a resize operations, since we have a new element
-    A[dim].resize(A[dim].size()+1);
-    
-    // update :)
-    for( int i = A[dim].size()-1; i >= idx; i--)
-    {
-      if (pf != NULL)
-        A[dim][i] = A[dim][i-1] + pow(pf->f ( 1, it->second ), 2);
-      else
-        A[dim][i] = A[dim][i-1] + pow(it->second, 2);
-    }   
-  }
-}
-
-void FastMinKernel::updateLookupTableForKVNApproximation(const NICE::SparseVector & _v, double * T, const Quantization & q, const ParameterizedFunction *pf) const
-{
-  if (T == NULL)
-  {
-    fthrow(Exception, "FastMinKernel::updateLookupTableForKernelVectorNorm LUT not initialized, run FastMinKernel::hikPrepareLookupTableForKernelVectorNorm first!");
-    return;
-  }
-  
-  // number of quantization bins
-  uint hmax = q.size();
-
-  // store (transformed) prototypes
-  double *prototypes = new double [ hmax ];
-  for ( uint i = 0 ; i < hmax ; i++ )
-    if ( pf != NULL ) {
-      // FIXME: the transformed prototypes could change from dimension to another dimension
-      // We skip this flexibility ...but it should be changed in the future
-      prototypes[i] = pf->f ( 1, q.getPrototype(i) );
-    } else {
-      prototypes[i] = q.getPrototype(i);
-    }
-   
-  // loop through all dimensions
-  for (NICE::SparseVector::const_iterator it = _v.begin(); it != _v.end(); it++)
-  {
-    int dim(it->first);
-
-    double x_i = it->second;
-    //as usually, we quantize the original features, but use the quantized transformed features lateron
-    int q_bin = q.quantize(x_i);      
-
-    //TODO we could speed up this with first do a binary search for the position where the min changes, and then do two separate for-loops
-    for (uint j = 0; j < hmax; j++)
-    {
-      double fval;
-      
-      if (q_bin > j)
-        fval = prototypes[j]; //the prototypes are already transformed
-      else
-      {
-        if (pf != NULL)
-          fval = pf->f( 1, x_i );
-        else
-          fval = x_i;
-      }
-      
-      // pay attention: we use either the quantized prototypes or the REAL feature values, not the quantized ones!
-      T[ dim*hmax + j ] += pow( fval, 2 );
-    }
-  }
-  
-  delete [] prototypes;  
-}

FastMinKernel.h

@@ -352,6 +352,10 @@ namespace NICE {
       */
       double* hikPrepareLookupTableForKVNApproximation(const Quantization & q, const ParameterizedFunction *pf = NULL) const;
       
+    //////////////////////////////////////////
+    // variance computation: sparse inputs
+    //////////////////////////////////////////      
+      
       /**
       * @brief Approximate norm = |k_*|^2 using the minimum kernel trick and exploiting sparsity of the given feature vector. Approximation does not considere mixed terms between dimensions.
       * @author Alexander Freytag
@@ -386,76 +390,50 @@ namespace NICE {
       */      
       void hikComputeKernelVector( const NICE::SparseVector & xstar, NICE::Vector & kstar) const;
       
-      /** Persistent interface */
-      virtual void restore ( std::istream & is, int format = 0 );
-      virtual void store ( std::ostream & os, int format = 0 ) const; 
-      virtual void clear ();
-      
-      // ----------------- INCREMENTAL LEARNING METHODS -----------------------
+    //////////////////////////////////////////
+    // variance computation: non-sparse inputs
+    //////////////////////////////////////////     
       
       /**
-      * @brief Add a new example to the feature-storage. You have to update the corresponding variables explicitely after that.
-      * @author Alexander Freytag
-      * @date 25-04-2012 (dd-mm-yyyy)
-      *
-      * @param _v new feature vector
-      */       
-      void addExample(const NICE::SparseVector & _v, const ParameterizedFunction *pf = NULL);
-      /**
-      * @brief Add a new example to the feature-storage. You have to update the corresponding variables explicitely after that.
-      * @author Alexander Freytag
-      * @date 25-04-2012 (dd-mm-yyyy)
-      *
-      * @param _v new feature vector
-      */       
-      void addExample(const std::vector<double> & _v, const ParameterizedFunction *pf = NULL);
-      
-      /**
-      * @brief Updates A and B matrices for fast kernel multiplications and kernel sums. You need to compute the new alpha value and run addExample first!
-      * @author Alexander Freytag
-      * @date 25-04-2012 (dd-mm-yyyy)
-      *
-      * @param _v new feature vector
-      * @param alpha new alpha value for the corresponding feature
-      * @param A precomputed matrix A which will be updated accordingly
-      * @param B precomputed matrix B which will be updated accordingly
-      * @param pf optional feature transformation
-      */       
-      void updatePreparationForAlphaMultiplications(const NICE::SparseVector & _v, const double & alpha, NICE::VVector & A, NICE::VVector & B, const ParameterizedFunction *pf = NULL) const;
-      /**
-      * @brief Updates LUT T for very fast kernel multiplications and kernel sums. You need to compute the new alpha value and run addExample first!
+      * @brief Approximate norm = |k_*|^2 using the minimum kernel trick and exploiting sparsity of the given feature vector. Approximation does not considere mixed terms between dimensions.
       * @author Alexander Freytag
-      * @date 26-04-2012 (dd-mm-yyyy)
-      *
-      * @param _v new feature vector
-      * @param alpha new alpha value for the corresponding feature
-      * @param T precomputed lookup table, which will be updated
-      * @param q quantization object to quantize possible test samples
+      * @date 19-12-2013 (dd-mm-yyyy)
+      * 
+      * @param A pre-computation matrix (VVector) (use the prepare method) 
+      * @param xstar new feature vector (Vector)
+      * @param norm result of the squared norm approximation
       * @param pf optional feature transformation
-      */       
-      void updateLookupTableForAlphaMultiplications(const NICE::SparseVector & _v, const double & alpha, double * T, const Quantization & q, const ParameterizedFunction *pf = NULL) const;
+      */
+      void hikComputeKVNApproximation(const NICE::VVector & A, const NICE::Vector & xstar, double & norm, const ParameterizedFunction *pf = NULL ) ;
       
       /**
-      * @brief Updates matrix A for approximations of the kernel vector norm. You need to run addExample first!
+      * @brief Approximate norm = |k_*|^2 using a large lookup table created by hikPrepareSquaredKernelVector and hikPrepareSquaredKernelVectorFast or directly using hikPrepareLookupTableForSquaredKernelVector. Approximation does not considere mixed terms between dimensions.
       * @author Alexander Freytag
-      * @date 26-04-2012 (dd-mm-yyyy)
+      * @date 19-12-2013 (dd-mm-yyyy)
       *
-      * @param _v new feature vector
-      * @param A precomputed matrix A which will be updated accordingly
-      * @param pf optional feature transformation
-      */       
-      void updatePreparationForKVNApproximation(const NICE::SparseVector & _v, NICE::VVector & A, const ParameterizedFunction *pf = NULL) const;
+      * @param Tlookup large lookup table
+      * @param q Quantization object
+      * @param xstar feature vector (indirect k_*)
+      * @param norm result of the calculation
+      */
+      void hikComputeKVNApproximationFast(const double *Tlookup, const Quantization & q, const NICE::Vector & xstar, double & norm ) const;      
+      
       /**
-      * @brief Updates LUT T for fast approximations of the kernel vector norm. You need to run addExample first!
+      * @brief Compute the kernel vector k_* between training examples and test example. Runtime. O(n \times D). Does not exploit sparsity - deprecated!
       * @author Alexander Freytag
-      * @date 26-04-2012 (dd-mm-yyyy)
+      * @date 19-12-2013 (dd-mm-yyyy)
       *
-      * @param _v new feature vector
-      * @param T precomputed lookup table, which will be updated
-      * @param q quantization object to quantize possible test samples
-      * @param pf optional feature transformation
-      */       
-      void updateLookupTableForKVNApproximation(const NICE::SparseVector & _v, double * T, const Quantization & q, const ParameterizedFunction *pf = NULL) const;
+      * @param xstar feature vector
+      * @param kstar kernel vector
+      */      
+      void hikComputeKernelVector( const NICE::Vector & xstar, NICE::Vector & kstar) const;      
+      
+      /** Persistent interface */
+      virtual void restore ( std::istream & is, int format = 0 );
+      virtual void store ( std::ostream & os, int format = 0 ) const; 
+      virtual void clear ();
+      
+     
 
   };
 

GMHIKernel.cpp

@@ -198,12 +198,3 @@ void GMHIKernel::setApproximationScheme(const int & _approxScheme)
 {
   this->fmk->setApproximationScheme(_approxScheme);
 }
-
-// ----------------- INCREMENTAL LEARNING METHODS -----------------------
-void GMHIKernel::addExample(const NICE::SparseVector & x, const NICE::Vector & binLabels)
-{
-  // we could add the example to the fmk, but we won't do it here
-  // reason: if we have a balanced learning, we have multiple identical GMHI-objects
-  // if we would add the example here, it would be added as often as we have those objects
-  // therefor we add the example already in the FMKGPHypOpt class
-}

GMHIKernel.h

@@ -78,8 +78,6 @@ class GMHIKernel : public ImplicitKernelMatrix
     virtual void store ( std::ostream & os, int format = 0 ) const {};//fmk->store( os, format );};
     virtual void clear () {};
     
-    virtual void addExample(const NICE::SparseVector & x, const NICE::Vector & binLabels);
-    
     void setFastMinKernel(NICE::FastMinKernel * _fmk){fmk = _fmk;};
      
 };

GPHIKClassifier.cpp

@@ -30,7 +30,7 @@ GPHIKClassifier::GPHIKClassifier( const Config *conf, const string & confSection
   
   if ( conf == NULL )
   {
-     fthrow(Exception, "GPHIKClassifier: the config is NULL -- use a default config and the restore-function instaed!");
+     fthrow(Exception, "GPHIKClassifier: the config is NULL -- use a default config or the restore-function instaed!");
   }
   else
     this->init(conf, confSection);
@@ -53,8 +53,6 @@ void GPHIKClassifier::init(const Config *conf, const string & confSection)
   double parameterLowerBound = conf->gD(confSection, "parameter_lower_bound", 1.0 );
   double parameterUpperBound = conf->gD(confSection, "parameter_upper_bound", 5.0 );
 
-  if (gphyper == NULL)
-    this->gphyper = new FMKGPHyperparameterOptimization;
   this->noise = conf->gD(confSection, "noise", 0.01);
 
   string transform = conf->gS(confSection, "transform", "absexp" );
@@ -90,16 +88,16 @@ void GPHIKClassifier::init(const Config *conf, const string & confSection)
   }
    
   //how do we approximate the predictive variance for classification uncertainty?
-  string varianceApproximationString = conf->gS(confSection, "varianceApproximation", "approximate_fine"); //default: fine approximative uncertainty prediction
-  if ( (varianceApproximationString.compare("approximate_rough") == 0) || ((varianceApproximationString.compare("1") == 0)) )
+  string s_varianceApproximation = conf->gS(confSection, "varianceApproximation", "approximate_fine"); //default: fine approximative uncertainty prediction
+  if ( (s_varianceApproximation.compare("approximate_rough") == 0) || ((s_varianceApproximation.compare("1") == 0)) )
   {
     this->varianceApproximation = APPROXIMATE_ROUGH;
   }
-  else if ( (varianceApproximationString.compare("approximate_fine") == 0) || ((varianceApproximationString.compare("2") == 0)) )
+  else if ( (s_varianceApproximation.compare("approximate_fine") == 0) || ((s_varianceApproximation.compare("2") == 0)) )
   {
     this->varianceApproximation = APPROXIMATE_FINE;
   }
-  else if ( (varianceApproximationString.compare("exact") == 0)  || ((varianceApproximationString.compare("3") == 0)) )
+  else if ( (s_varianceApproximation.compare("exact") == 0)  || ((s_varianceApproximation.compare("3") == 0)) )
   {
     this->varianceApproximation = EXACT;
   }
@@ -107,7 +105,7 @@ void GPHIKClassifier::init(const Config *conf, const string & confSection)
   {
     this->varianceApproximation = NONE;
   } 
-  std::cerr << "varianceApproximationStrategy: " << varianceApproximationString  << std::endl;
+  std::cerr << "varianceApproximationStrategy: " << s_varianceApproximation  << std::endl;
 }
 
 void GPHIKClassifier::classify ( const SparseVector * example,  int & result, SparseVector & scores )
@@ -124,6 +122,9 @@ void GPHIKClassifier::classify ( const NICE::Vector * example,  int & result, Sp
 
 void GPHIKClassifier::classify ( const SparseVector * example,  int & result, SparseVector & scores, double & uncertainty )
 {
+  if (gphyper == NULL)
+     fthrow(Exception, "Classifier not trained yet -- aborting!" );
+  
   scores.clear();
   
   int classno = gphyper->classify ( *example, scores );
@@ -138,9 +139,7 @@ void GPHIKClassifier::classify ( const SparseVector * example,  int & result, Sp
   {
     if (varianceApproximation != NONE)
     {
-      NICE::Vector uncertainties;
-      this->predictUncertainty( example, uncertainties );
-      uncertainty = uncertainties.Max();
+      this->predictUncertainty( example, uncertainty );
     }  
     else
     {
@@ -157,6 +156,9 @@ void GPHIKClassifier::classify ( const SparseVector * example,  int & result, Sp
 
 void GPHIKClassifier::classify ( const NICE::Vector * example,  int & result, SparseVector & scores, double & uncertainty )
 {
+  if (gphyper == NULL)
+     fthrow(Exception, "Classifier not trained yet -- aborting!" );  
+  
   scores.clear();
   
   int classno = gphyper->classify ( *example, scores );
@@ -166,13 +168,12 @@ void GPHIKClassifier::classify ( const NICE::Vector * example,  int & result, Sp
   }
   
   result = scores.maxElement();
-   
+  
   if (uncertaintyPredictionForClassification)
   {
     if (varianceApproximation != NONE)
     {
-      std::cerr << "ERROR: Uncertainty computation is currently not supported for NICE::Vector - use SparseVector instead" << std::endl;
-      uncertainty = std::numeric_limits<double>::max();
+      this->predictUncertainty( example, uncertainty );
     }  
     else
     {
@@ -184,7 +185,7 @@ void GPHIKClassifier::classify ( const NICE::Vector * example,  int & result, Sp
   {
     //do nothing
     uncertainty = std::numeric_limits<double>::max();
-  }    
+  }  
 }
 
 /** training process */
@@ -201,6 +202,8 @@ void GPHIKClassifier::train ( const std::vector< NICE::SparseVector *> & example
   if (verbose)
     std::cerr << "Time used for setting up the fmk object: " << t.getLast() << std::endl;  
   
+  if (gphyper != NULL)
+     delete gphyper;
   gphyper = new FMKGPHyperparameterOptimization ( confCopy, pf, fmk, confSection ); 
 
   if (verbose)
@@ -209,14 +212,36 @@ void GPHIKClassifier::train ( const std::vector< NICE::SparseVector *> & example
   // go go go
   gphyper->optimize ( labels );
   if (verbose)
-    std::cerr << "optimization done, now prepare for the uncertainty prediction" << std::endl;
+    std::cerr << "optimization done" << std::endl;
   
-  if ( (varianceApproximation == APPROXIMATE_ROUGH) )
+  if ( ( varianceApproximation != NONE ) )
   {
-    //prepare for variance computation (approximative)
-    gphyper->prepareVarianceApproximation();
+    std::cerr << "now prepare for the uncertainty prediction" << std::endl;
+    
+    switch (varianceApproximation)    
+    {
+      case APPROXIMATE_ROUGH:
+      {
+        gphyper->prepareVarianceApproximationRough();
+        break;
+      }
+      case APPROXIMATE_FINE:
+      {
+        gphyper->prepareVarianceApproximationFine();
+        break;
+      }    
+      case EXACT:
+      {
+       //nothing to prepare
+        break;
+      }
+      default:
+      {
+       //nothing to prepare
+      }
+    }
   }
-  //for exact variance computation, we do not have to prepare anything
+
 
   // clean up all examples ??
   if (verbose)
@@ -236,6 +261,8 @@ void GPHIKClassifier::train ( const std::vector< SparseVector *> & examples, std
   if (verbose)
     std::cerr << "Time used for setting up the fmk object: " << t.getLast() << std::endl;  
   
+  if (gphyper != NULL)
+     delete gphyper;
   gphyper = new FMKGPHyperparameterOptimization ( confCopy, pf, fmk, confSection ); 
 
   if (verbose)
@@ -245,12 +272,33 @@ void GPHIKClassifier::train ( const std::vector< SparseVector *> & examples, std
   if (verbose)
     std::cerr << "optimization done, now prepare for the uncertainty prediction" << std::endl;
   
-  if ( (varianceApproximation == APPROXIMATE_ROUGH) )
+  if ( ( varianceApproximation != NONE ) )
   {
-    //prepare for variance computation (approximative)
-    gphyper->prepareVarianceApproximation();
+    std::cerr << "now prepare for the uncertainty prediction" << std::endl;
+    
+    switch (varianceApproximation)    
+    {
+      case APPROXIMATE_ROUGH:
+      {
+        gphyper->prepareVarianceApproximationRough();
+        break;
+      }
+      case APPROXIMATE_FINE:
+      {
+        gphyper->prepareVarianceApproximationFine();
+        break;
+      }    
+      case EXACT:
+      {
+       //nothing to prepare
+        break;
+      }
+      default:
+      {
+       //nothing to prepare
+      }
+    }
   }
-  //for exact variance computation, we do not have to prepare anything
 
   // clean up all examples ??
   if (verbose)
@@ -271,32 +319,67 @@ GPHIKClassifier *GPHIKClassifier::clone () const
   return NULL;
 }
   
-void GPHIKClassifier::predictUncertainty( const NICE::SparseVector * example, NICE::Vector & uncertainties )
+void GPHIKClassifier::predictUncertainty( const NICE::SparseVector * example, double & uncertainty )
 {  
+  if (gphyper == NULL)
+     fthrow(Exception, "Classifier not trained yet -- aborting!" );  
+  
   //we directly store the predictive variances in the vector, that contains the classification uncertainties lateron to save storage
   switch (varianceApproximation)    
   {
     case APPROXIMATE_ROUGH:
     {
-      gphyper->computePredictiveVarianceApproximateRough( *example, uncertainties );
+      gphyper->computePredictiveVarianceApproximateRough( *example, uncertainty );
       break;
     }
     case APPROXIMATE_FINE:
     {
-      gphyper->computePredictiveVarianceApproximateFine( *example, uncertainties );
+        std::cerr << "predict uncertainty fine" << std::endl;
+      gphyper->computePredictiveVarianceApproximateFine( *example, uncertainty );
       break;
     }    
     case EXACT:
     {
-      gphyper->computePredictiveVarianceExact( *example, uncertainties );
+      gphyper->computePredictiveVarianceExact( *example, uncertainty );
       break;
     }
     default:
     {
-//       std::cerr << "No Uncertainty Prediction at all" << std::endl;
       fthrow(Exception, "GPHIKClassifier - your settings disabled the variance approximation needed for uncertainty prediction.");
-//       uncertainties.resize( 1 );
-//       uncertainties.set( numeric_limits<double>::max() );
+//       uncertainty = numeric_limits<double>::max();
+//       break;
+    }
+  }
+}
+
+void GPHIKClassifier::predictUncertainty( const NICE::Vector * example, double & uncertainty )
+{  
+  if (gphyper == NULL)
+     fthrow(Exception, "Classifier not trained yet -- aborting!" );  
+  
+  //we directly store the predictive variances in the vector, that contains the classification uncertainties lateron to save storage
+  switch (varianceApproximation)    
+  {
+    case APPROXIMATE_ROUGH:
+    {
+      gphyper->computePredictiveVarianceApproximateRough( *example, uncertainty );
+      break;
+    }
+    case APPROXIMATE_FINE:
+    {
+        std::cerr << "predict uncertainty fine" << std::endl;
+      gphyper->computePredictiveVarianceApproximateFine( *example, uncertainty );
+      break;
+    }    
+    case EXACT:
+    {
+      gphyper->computePredictiveVarianceExact( *example, uncertainty );
+      break;
+    }
+    default:
+    {
+      fthrow(Exception, "GPHIKClassifier - your settings disabled the variance approximation needed for uncertainty prediction.");
+//       uncertainty = numeric_limits<double>::max();
 //       break;
     }
   }
@@ -355,6 +438,9 @@ void GPHIKClassifier::restore ( std::istream & is, int format )
 
 void GPHIKClassifier::store ( std::ostream & os, int format ) const
 {
+  if (gphyper == NULL)
+     fthrow(Exception, "Classifier not trained yet -- aborting!" );  
+  
   if (os.good())
   {
     os.precision (numeric_limits<double>::digits10 + 1);
@@ -379,16 +465,10 @@ void GPHIKClassifier::store ( std::ostream & os, int format ) const
   }
 }
 
-void GPHIKClassifier::addExample( const NICE::SparseVector * example, const double & label, const bool & performOptimizationAfterIncrement)
-{
-  gphyper->addExample( *example, label, performOptimizationAfterIncrement );
-}
-
-void GPHIKClassifier::addMultipleExamples( const std::vector< const NICE::SparseVector *> & newExamples, const NICE::Vector & newLabels, const bool & performOptimizationAfterIncrement)
+std::set<int> GPHIKClassifier::getKnownClassNumbers ( ) const
 {
-  //are new examples available? If not, nothing has to be done
-  if ( newExamples.size() < 1)
-    return;
+  if (gphyper == NULL)
+     fthrow(Exception, "Classifier not trained yet -- aborting!" );  
   
-  gphyper->addMultipleExamples( newExamples, newLabels, performOptimizationAfterIncrement );
+  return gphyper->getKnownClassNumbers();
 }

GPHIKClassifier.h

@@ -109,7 +109,7 @@ class GPHIKClassifier
      * @param example (non-sparse Vector) to be classified given in a non-sparse representation
      * @param result (int) class number of most likely class
      * @param scores (SparseVector) classification scores for known classes
-     * @param uncertainty (double*) predictive variance of the classification result, if computed
+     * @param uncertainty (double) predictive variance of the classification result, if computed
      */    
     void classify ( const NICE::Vector * example,  int & result, NICE::SparseVector & scores, double & uncertainty );    
 
@@ -143,12 +143,21 @@ class GPHIKClassifier
      * @date 19-06-2012 (dd-mm-yyyy)
      * @author Alexander Freytag
      * @param examples example for which the classification uncertainty shall be predicted, given in a sparse representation
-     * @param uncertainties contains the resulting classification uncertainties (1 entry for standard setting, m entries for binary-balanced setting)
+     * @param uncertainty contains the resulting classification uncertainty
      */       
-    void predictUncertainty( const NICE::SparseVector * example, NICE::Vector & uncertainties );
+    void predictUncertainty( const NICE::SparseVector * example, double & uncertainty );
     
-    void addExample( const NICE::SparseVector * example, const double & label, const bool & performOptimizationAfterIncrement = true);
-    void addMultipleExamples( const std::vector< const NICE::SparseVector * > & newExamples, const NICE::Vector & newLabels, const bool & performOptimizationAfterIncrement = true);
+    /** 
+     * @brief prediction of classification uncertainty
+     * @date 19-12-2013 (dd-mm-yyyy)
+     * @author Alexander Freytag
+     * @param examples example for which the classification uncertainty shall be predicted, given in a non-sparse representation
+     * @param uncertainty contains the resulting classification uncertainty
+     */       
+    void predictUncertainty( const NICE::Vector * example, double & uncertainty );    
+    
+    std::set<int> getKnownClassNumbers ( ) const;
+
 };
 
 }

GPLikelihoodApprox.cpp

@@ -52,21 +52,13 @@ GPLikelihoodApprox::GPLikelihoodApprox( const map<int, Vector> & binaryLabels,
   this->verifyApproximation = verifyApproximation;
   
   this->nrOfEigenvaluesToConsider = _nrOfEigenvaluesToConsider;
-  
-  lastAlphas = NULL;
-  
+    
   this->verbose = false;
   this->debug = false;
-  
-  this->usePreviousAlphas = true;
-
 }
 
 GPLikelihoodApprox::~GPLikelihoodApprox()
 {
-  //delete the pointer, but not the content (which is stored somewhere else)
-  if (lastAlphas != NULL)
-    lastAlphas = NULL;  
 }
 
 void GPLikelihoodApprox::calculateLikelihood ( double mypara, const FeatureMatrix & f, const std::map< int, NICE::Vector > & yset, double noise, double lambdaMax )
@@ -98,9 +90,9 @@ void GPLikelihoodApprox::calculateLikelihood ( double mypara, const FeatureMatri
   cerr << "chol * chol^T: " << ( choleskyMatrix * choleskyMatrix.transpose() )(0,0,4,4) << endl;
 
   double gt_dataterm = 0.0;
-  for ( map< int, NICE::Vector >::const_iterator i = yset.begin(); i != yset.end(); i++ )
+  for ( std::map< int, NICE::Vector >::const_iterator i = yset.begin(); i != yset.end(); i++ )
   {
-    const Vector & y = i->second;
+    const NICE::Vector & y = i->second;
     Vector gt_alpha;
     choleskySolve ( choleskyMatrix, y, gt_alpha );
     cerr << "cholesky error: " << (K*gt_alpha - y).normL2() << endl;
@@ -113,31 +105,23 @@ void GPLikelihoodApprox::calculateLikelihood ( double mypara, const FeatureMatri
   cerr << "Something of K: " << K(0, 0, 4, 4) << endl;
   cerr << "frob norm: gt:" << K.frobeniusNorm() << endl;
   
-  /*try {
-    Vector *eigenv = eigenvalues ( K ); 
-    cerr << "lambda_max: gt:" << eigenv->Max() << " est:" << lambdaMax << endl; 
-    delete eigenv;
-  } catch (...) {
-    cerr << "NICE eigenvalues function failed!" << endl;
-  }*/
-
+  
   double gt_nlikelihood = gt_logdet + gt_dataterm;
   cerr << "OPTGT: " << mypara << " " << gt_nlikelihood << " " << gt_logdet << " " << gt_dataterm << endl;
 }
 
-void GPLikelihoodApprox::computeAlphaDirect(const OPTIMIZATION::matrix_type & x)
+void GPLikelihoodApprox::computeAlphaDirect(const OPTIMIZATION::matrix_type & x, const NICE::Vector & eigenValues )
 {
   Timer t;
-//   NICE::Vector diagonalElements;
   
-//   ikm->getDiagonalElements ( diagonalElements );
+  NICE::Vector diagonalElements; 
+  ikm->getDiagonalElements ( diagonalElements );
 
   // set simple jacobi pre-conditioning
   ILSConjugateGradients *linsolver_cg = dynamic_cast<ILSConjugateGradients *> ( linsolver );
 
-//   //TODO why do we need this?  
-//   if ( linsolver_cg != NULL )
-//     linsolver_cg->setJacobiPreconditioner ( diagonalElements );
+  if ( linsolver_cg != NULL )
+    linsolver_cg->setJacobiPreconditioner ( diagonalElements );
   
 
   // all alpha vectors will be stored!
@@ -153,8 +137,6 @@ void GPLikelihoodApprox::computeAlphaDirect(const OPTIMIZATION::matrix_type & x)
     if (verbose)
     {
       std::cerr << "Solving linear equation system for class " << classCnt << " ..." << std::endl;
-      std::cerr << "Size of the kernel matrix " << ikm->rows() << std::endl;
-      std::cerr << "binary label: " << j->second << std::endl;
     }
 
     /** About finding a good initial solution
@@ -171,23 +153,10 @@ void GPLikelihoodApprox::computeAlphaDirect(const OPTIMIZATION::matrix_type & x)
      * v = y ....which is somehow a weird assumption (cf Kernel PCA)
      *  This reduces the number of iterations by 5 or 8
      */
-    Vector alpha;
+    NICE::Vector alpha;
     
-    if ( (usePreviousAlphas) && (lastAlphas != NULL) )
-    {
-      std::map<int, NICE::Vector>::iterator alphaIt = lastAlphas->begin();
-      alpha = (*lastAlphas)[classCnt];
-    }
-    else  
-    {
-      //TODO hand over the eigenmax
-      alpha = (binaryLabels[classCnt] ); //* (1.0 / eigenmax[0]) );
-    }
+    alpha = (binaryLabels[classCnt] * (1.0 / eigenValues[0]) );
     
-    NICE::Vector initialAlpha;
-    if ( verbose )
-     initialAlpha = alpha;
-
     if ( verbose )
       std::cerr << "Using the standard solver ..." << std::endl;
 
@@ -205,8 +174,8 @@ void GPLikelihoodApprox::computeAlphaDirect(const OPTIMIZATION::matrix_type & x)
 
 double GPLikelihoodApprox::evaluate(const OPTIMIZATION::matrix_type & x)
 {
-  Vector xv;
- 
+  NICE::Vector xv;
+   
   xv.resize ( x.rows() );
   for ( uint i = 0 ; i < x.rows(); i++ )
     xv[i] = x(i,0);
@@ -215,7 +184,7 @@ double GPLikelihoodApprox::evaluate(const OPTIMIZATION::matrix_type & x)
   unsigned long hashValue = xv.getHashValue();
   if (verbose)  
     std::cerr << "Current parameter: " << xv << " (weird hash value is " << hashValue << ")" << std::endl;
-  map<unsigned long, double>::const_iterator k = alreadyVisited.find(hashValue);
+  std::map<unsigned long, double>::const_iterator k = alreadyVisited.find(hashValue);
   
   if ( k != alreadyVisited.end() )
   {
@@ -244,8 +213,8 @@ double GPLikelihoodApprox::evaluate(const OPTIMIZATION::matrix_type & x)
   if (verbose)  
     std::cerr << "Calculating eigendecomposition " << ikm->rows() << " x " << ikm->cols() << std::endl;
   t.start();
-  Vector eigenmax;
-  Matrix eigenmaxvectors;
+  NICE::Vector eigenmax;
+  NICE::Matrix eigenmaxvectors;
  
   int rank = nrOfEigenvaluesToConsider;
 
@@ -255,9 +224,7 @@ double GPLikelihoodApprox::evaluate(const OPTIMIZATION::matrix_type & x)
   // the current implementation converges very quickly
   //old version: just use the first eigenvalue
   
-  //NOTE
-  // in theory, we have these values already on hand since we've done it in FMKGPHypOpt.
-  // Think about wether to give them as input to this function or not
+  // we have to re-compute EV and EW in all cases, since we change the hyper parameter and thereby the kernel matrix 
   eig->getEigenvalues( *ikm, eigenmax, eigenmaxvectors, rank ); 
   if (verbose)
     std::cerr << "eigenmax: " << eigenmax << std::endl;
@@ -278,12 +245,12 @@ double GPLikelihoodApprox::evaluate(const OPTIMIZATION::matrix_type & x)
   
 
   // all alpha vectors will be stored!
-  map<int, Vector> alphas;
+  std::map<int, NICE::Vector> alphas;
 
   // This has to be done m times for the multi-class case
   if (verbose)
     std::cerr << "run ILS for every bin label. binaryLabels.size(): " << binaryLabels.size() << std::endl;
-  for ( map<int, Vector>::const_iterator j = binaryLabels.begin(); j != binaryLabels.end() ; j++)
+  for ( std::map<int, NICE::Vector>::const_iterator j = binaryLabels.begin(); j != binaryLabels.end() ; j++)
   {
     // (b) y^T (K+sI)^{-1} y
     int classCnt = j->first;
@@ -307,22 +274,10 @@ double GPLikelihoodApprox::evaluate(const OPTIMIZATION::matrix_type & x)
      * v = y ....which is somehow a weird assumption (cf Kernel PCA)
      *  This reduces the number of iterations by 5 or 8
      */
-    Vector alpha;
+    NICE::Vector alpha;
     
-    if ( (usePreviousAlphas) && (lastAlphas != NULL) )
-    {
-      std::map<int, NICE::Vector>::iterator alphaIt = lastAlphas->begin();
-      alpha = (*lastAlphas)[classCnt];
-    }
-    else  
-    {
-      alpha = (binaryLabels[classCnt] * (1.0 / eigenmax[0]) );
-    }
+    alpha = (binaryLabels[classCnt] * (1.0 / eigenmax[0]) );
     
-    Vector initialAlpha;
-    if ( verbose )
-     initialAlpha = alpha;
-
     if ( verbose )
       cerr << "Using the standard solver ..." << endl;
 
@@ -330,22 +285,6 @@ double GPLikelihoodApprox::evaluate(const OPTIMIZATION::matrix_type & x)
     linsolver->solveLin ( *ikm, binaryLabels[classCnt], alpha );
     t.stop();
    
-    //TODO This is only important for the incremental learning stuff.
-//     if ( verbose )
-//     {
-//       double initialAlphaNorm ( initialAlpha.normL1() );
-//       //compute the difference
-//       initialAlpha -= alpha;
-//       //take the abs of the differences
-//       initialAlpha.absInplace();
-//       //and compute a final score using a suitable norm
-// //       double difference( initialAlpha.normInf() );
-//       double difference( initialAlpha.normL1() );
-//       std::cerr << "debug -- last entry of new alpha: " << abs(alpha[alpha.size() -1 ]) << std::endl;
-//       std::cerr << "debug -- difference using inf norm: " << difference  << std::endl;
-//       std::cerr << "debug -- relative difference using inf norm: " << difference / initialAlphaNorm  << std::endl;
-//     }
-
 
     if ( verbose )
       std::cerr << "Time used for solving (K + sigma^2 I)^{-1} y: " << t.getLast() << std::endl;
@@ -418,21 +357,12 @@ void GPLikelihoodApprox::setParameterUpperBound(const double & _parameterUpperBo
   parameterUpperBound = _parameterUpperBound;
 }
 
-void GPLikelihoodApprox::setLastAlphas(std::map<int, NICE::Vector> * _lastAlphas)
-{
-  lastAlphas = _lastAlphas;
-}
 
 void GPLikelihoodApprox::setBinaryLabels(const std::map<int, Vector> & _binaryLabels)
 {
   binaryLabels = _binaryLabels;
 }
 
-void GPLikelihoodApprox::setUsePreviousAlphas( const bool & _usePreviousAlphas )
-{
-  this->usePreviousAlphas = _usePreviousAlphas; 
-}
-
 void GPLikelihoodApprox::setVerbose( const bool & _verbose )
 {
   this->verbose = _verbose;

GPLikelihoodApprox.h

@@ -111,7 +111,7 @@ class GPLikelihoodApprox : public OPTIMIZATION::CostFunction
     *
     * @return void
     */    
-    void computeAlphaDirect(const OPTIMIZATION::matrix_type & x);
+    void computeAlphaDirect(const OPTIMIZATION::matrix_type & x, const NICE::Vector & eigenValues);
     
     /**
     * @brief Evaluate the likelihood for given hyperparameters
@@ -124,18 +124,19 @@ class GPLikelihoodApprox : public OPTIMIZATION::CostFunction
      
     
     // ------ get and set methods ------
-    const Vector & getBestParameters () const { return min_parameter; };
+    const NICE::Vector & getBestParameters () const { return min_parameter; };
     const std::map<int, Vector> & getBestAlphas () const { return min_alphas; };
     
     void setParameterLowerBound(const double & _parameterLowerBound);
     void setParameterUpperBound(const double & _parameterUpperBound);
     
-    void setLastAlphas(std::map<int, NICE::Vector> * _lastAlphas);
     void setBinaryLabels(const std::map<int, Vector> & _binaryLabels);
     
-    void setUsePreviousAlphas( const bool & _usePreviousAlphas );
     void setVerbose( const bool & _verbose );
     void setDebug( const bool & _debug );
+    
+    bool getVerbose ( ) { return verbose; } ;
+    bool getDebug ( ) { return debug; } ;
 };
 
 }

IKMLinearCombination.cpp

@@ -218,13 +218,3 @@ ImplicitKernelMatrix * IKMLinearCombination::getModel(const uint & idx) const
   else
     return NULL;
 }
-
-// ----------------- INCREMENTAL LEARNING METHODS -----------------------
-void IKMLinearCombination::addExample(const NICE::SparseVector & x, const NICE::Vector & binLabels)
-{
-  for ( vector<ImplicitKernelMatrix *>::iterator i = matrices.begin(); i != matrices.end(); i++ )
-  {
-    ImplicitKernelMatrix *ikm = *i;
-    ikm->addExample(x, binLabels);
-  }
-}

IKMLinearCombination.h

@@ -71,8 +71,6 @@ class IKMLinearCombination : public ImplicitKernelMatrix
     virtual void restore ( std::istream & is, int format = 0 ) {};
     virtual void store ( std::ostream & os, int format = 0 ) const {};  
     virtual void clear () {};
-    
-    void addExample(const NICE::SparseVector & x, const NICE::Vector & binLabels);
 
 };
 

IKMNoise.cpp

@@ -230,17 +230,3 @@ void IKMNoise::store ( std::ostream & os, int format ) const
   os << "nn: " << nn << std::endl;
   os << "labels: " << labels << std::endl;
 }
-
-// ----------------- INCREMENTAL LEARNING METHODS -----------------------
-void IKMNoise::addExample(const NICE::SparseVector & x, const NICE::Vector & binLabels)
-{
-  ++size;
-  if ( (np != 0) && (nn != 0) )
-  {
-    labels = binLabels;
-    if (binLabels[binLabels.size()-1] == 1)
-      ++np;
-    else
-      ++nn;
-  }
-}

IKMNoise.h

@@ -75,7 +75,6 @@ class IKMNoise : public ImplicitKernelMatrix
     virtual void store ( std::ostream & os, int format = 0 ) const; 
     virtual void clear () {};
     
-    void addExample(const NICE::SparseVector & x, const NICE::Vector & binLabels);
 
 };
 

ImplicitKernelMatrix.h

@@ -57,10 +57,9 @@ class ImplicitKernelMatrix : public GenericMatrix, NICE::Persistent
     virtual void clear () = 0;
     
     //high order methods
-    virtual void addExample(const NICE::SparseVector & x, const NICE::Vector & binLabels) = 0;
     virtual void  multiply (NICE::Vector &y, const NICE::Vector &x) const = 0;
 };
 
 }
 
-#endif
+#endif

matlab/GPHIK.cpp

@@ -0,0 +1,848 @@
+#include <math.h>
+#include <matrix.h>
+#include "mex.h"
+#include "classHandleMtoC.h"
+
+// NICE-core includes
+#include <core/basics/Config.h>
+#include <core/basics/Timer.h>
+#include <core/vector/MatrixT.h>
+#include <core/vector/VectorT.h>
+
+// gp-hik-core includes
+#include "gp-hik-core/GPHIKClassifier.h"
+
+using namespace std; //C basics
+using namespace NICE;  // nice-core
+
+/* Pass analyze_sparse a pointer to a sparse mxArray.  A sparse mxArray
+   only stores its nonzero elements.  The values of the nonzero elements 
+   are stored in the pr and pi arrays.  The tricky part of analyzing
+   sparse mxArray's is figuring out the indices where the nonzero
+   elements are stored.  (See the mxSetIr and mxSetJc reference pages
+   for details. */  
+std::vector< NICE::SparseVector * > convertSparseMatrixToNice(const mxArray *array_ptr)
+{
+  double  *pr;//, *pi;
+  mwIndex  *ir, *jc;
+  mwSize      col, total=0;
+  mwIndex   starting_row_index, stopping_row_index, current_row_index;
+  mwSize      i_numExamples, i_numDim;
+  
+  /* Get the starting positions of all four data arrays. */ 
+  pr = mxGetPr(array_ptr);
+//   pi = mxGetPi(array_ptr);
+  ir = mxGetIr(array_ptr);
+  jc = mxGetJc(array_ptr);
+  
+  // dimenions of the matrix -> feature dimension and number of examples
+  i_numExamples = mxGetM(array_ptr);  
+  i_numDim = mxGetN(array_ptr);
+    
+  // initialize output variable
+  std::vector< NICE::SparseVector * > sparseMatrix;
+  sparseMatrix.resize ( i_numExamples );
+    
+  for ( std::vector< NICE::SparseVector * >::iterator matIt = sparseMatrix.begin(); 
+        matIt != sparseMatrix.end(); matIt++)
+  {
+      *matIt = new NICE::SparseVector( i_numDim );
+  }  
+  
+  // now copy the data
+  for (col=0; col < i_numDim; col++)
+  { 
+    starting_row_index = jc[col]; 
+    stopping_row_index = jc[col+1]; 
+    
+    // empty column?
+    if (starting_row_index == stopping_row_index)
+      continue;
+    else
+    {
+      for ( current_row_index = starting_row_index; 
+            current_row_index < stopping_row_index; 
+	        current_row_index++)
+      {
+          //note: no complex data supported her
+          sparseMatrix[ ir[current_row_index] ]->insert( std::pair<int, double>( col, pr[total++] ) );
+      } // for-loop
+      
+    }
+  } // for-loop over columns
+  
+  return sparseMatrix;
+}
+
+
+// b_adaptIndexMtoC: if true, dim k will be inserted as k, not as k-1 (which would be the default for  M->C)
+NICE::SparseVector convertSparseVectorToNice(const mxArray* array_ptr, const bool & b_adaptIndexMtoC = false )
+{
+  double  *pr, *pi;
+  mwIndex  *ir, *jc;
+  mwSize      col, total=0;
+  mwIndex   starting_row_index, stopping_row_index, current_row_index;
+  mwSize      dimy, dimx;
+  
+  /* Get the starting positions of all four data arrays. */ 
+  pr = mxGetPr(array_ptr);
+  pi = mxGetPi(array_ptr);
+  ir = mxGetIr(array_ptr);
+  jc = mxGetJc(array_ptr);
+  
+  // dimenions of the matrix -> feature dimension and number of examples
+  dimy = mxGetM(array_ptr);  
+  dimx = mxGetN(array_ptr);
+  
+  double* ptr = mxGetPr(array_ptr);
+
+  if(dimx != 1 && dimy != 1)
+    mexErrMsgIdAndTxt("mexnice:error","Vector expected");
+  
+
+  NICE::SparseVector svec( std::max(dimx, dimy) );
+   
+  
+  if ( dimx > 1)
+  {
+    for ( mwSize row=0; row < dimx; row++)
+    { 
+        // empty column?
+        if (jc[row] == jc[row+1])
+        {
+          continue;
+        }
+        else
+        {
+          //note: no complex data supported her
+            double value ( pr[total++] );
+            if ( b_adaptIndexMtoC ) 
+                svec.insert( std::pair<int, double>( row+1,  value ) );
+            else
+                svec.insert( std::pair<int, double>( row,  value ) );
+        }
+    } // for loop over cols      
+  }
+  else
+  {
+    mwSize numNonZero = jc[1]-jc[0];
+    
+    for ( mwSize colNonZero=0; colNonZero < numNonZero; colNonZero++)
+    {
+        //note: no complex data supported her
+        double value ( pr[total++] );
+        if ( b_adaptIndexMtoC ) 
+            svec.insert( std::pair<int, double>( ir[colNonZero]+1, value  ) );
+        else
+            svec.insert( std::pair<int, double>( ir[colNonZero], value  ) );
+    }          
+  }
+
+  return svec;
+}
+
+// b_adaptIndexCtoM: if true, dim k will be inserted as k, not as k+1 (which would be the default for C->M)
+mxArray* convertSparseVectorFromNice( const NICE::SparseVector & scores, const bool & b_adaptIndexCtoM = false)
+{
+    mxArray * matlabSparseVec = mxCreateSparse( scores.getDim() /*m*/, 1/*n*/, scores.size()/*nzmax*/, mxREAL);
+    
+    // To make the returned sparse mxArray useful, you must initialize the pr, ir, jc, and (if it exists) pi arrays.    
+    // mxCreateSparse allocates space for:
+    // 
+    // A pr array of length nzmax.
+    // A pi array of length nzmax, but only if ComplexFlag is mxCOMPLEX in C (1 in Fortran).
+    // An ir array of length nzmax.
+    // A jc array of length n+1.  
+  
+    double* prPtr = mxGetPr(matlabSparseVec);
+    mwIndex * ir = mxGetIr( matlabSparseVec );
+    
+    mwIndex * jc = mxGetJc( matlabSparseVec );
+    jc[1] = scores.size(); jc[0] = 0; 
+    
+    
+    mwSize cnt = 0;
+        
+    for ( NICE::SparseVector::const_iterator myIt = scores.begin(); myIt != scores.end(); myIt++, cnt++ )
+    {
+        // set index
+        if ( b_adaptIndexCtoM ) 
+            ir[cnt] = myIt->first-1;
+        else
+            ir[cnt] = myIt->first;
+        
+        // set value
+        prPtr[cnt] = myIt->second;
+    }
+    
+    return matlabSparseVec;
+}
+
+
+mxArray* convertMatrixFromNice(NICE::Matrix & niceMatrix)
+{
+	mxArray *matlabMatrix = mxCreateDoubleMatrix(niceMatrix.rows(),niceMatrix.cols(),mxREAL);
+	double* matlabMatrixPtr = mxGetPr(matlabMatrix);
+
+	for(int i=0; i<niceMatrix.rows(); i++)
+    {
+		for(int j=0; j<niceMatrix.cols(); j++)
+		{
+			matlabMatrixPtr[i + j*niceMatrix.rows()] = niceMatrix(i,j);
+		}
+    }
+	return matlabMatrix;
+}
+
+NICE::Matrix convertMatrixToNice(const mxArray* matlabMatrix)
+{
+	//todo: do not assume double
+
+  const mwSize *dims;
+  int dimx, dimy, numdims;
+    //figure out dimensions
+  dims = mxGetDimensions(matlabMatrix);
+  numdims = mxGetNumberOfDimensions(matlabMatrix);
+  dimy = (int)dims[0]; dimx = (int)dims[1];
+  double* ptr = mxGetPr(matlabMatrix);
+
+  NICE::Matrix niceMatrix(ptr, dimy, dimx, NICE::Matrix::external); 
+
+  return niceMatrix;
+}
+
+mxArray* convertVectorFromNice(NICE::Vector & niceVector)
+{
+	//cout << "start convertVectorFromNice" << endl;
+	mxArray *matlabVector = mxCreateDoubleMatrix(niceVector.size(), 1, mxREAL);
+	double* matlabVectorPtr = mxGetPr(matlabVector);
+
+	for(int i=0;i<niceVector.size(); i++)
+    {
+        matlabVectorPtr[i] = niceVector[i];
+    }
+	return matlabVector;
+}
+
+NICE::Vector convertVectorToNice(const mxArray* matlabMatrix)
+{
+	//todo: do not assume double
+
+  const mwSize *dims;
+  int dimx, dimy, numdims;
+    //figure out dimensions
+  dims = mxGetDimensions(matlabMatrix);
+  numdims = mxGetNumberOfDimensions(matlabMatrix);
+  dimy = (int)dims[0]; dimx = (int)dims[1];
+  double* ptr = mxGetPr(matlabMatrix);
+
+  if(dimx != 1 && dimy != 1)
+    mexErrMsgIdAndTxt("mexnice:error","Vector expected");
+
+  int dim = max(dimx, dimy);    
+
+  NICE::Vector niceVector(dim, 0.0);
+  
+  for(int i=0;i<dim;i++)
+  {
+      niceVector(i) = ptr[i];
+  }
+
+  return niceVector;
+}
+
+
+
+std::string convertMatlabToString(const mxArray *matlabString)
+{
+  if(!mxIsChar(matlabString))
+    mexErrMsgIdAndTxt("mexnice:error","Expected string");
+
+  char *cstring = mxArrayToString(matlabString);
+  std::string s(cstring);
+  mxFree(cstring);
+  return s;
+}
+
+
+int convertMatlabToInt32(const mxArray *matlabInt32)
+{
+  if(!mxIsInt32(matlabInt32))
+    mexErrMsgIdAndTxt("mexnice:error","Expected int32");
+
+  int* ptr = (int*)mxGetData(matlabInt32);
+  return ptr[0];
+}
+
+double convertMatlabToDouble(const mxArray *matlabDouble)
+{
+  if(!mxIsDouble(matlabDouble))
+    mexErrMsgIdAndTxt("mexnice:error","Expected double");
+
+  double* ptr = (double*)mxGetData(matlabDouble);
+  return ptr[0];
+}
+
+Config parseParameters(const mxArray *prhs[], int nrhs)
+{
+  Config conf;
+  for(int i=0;i<nrhs;i+=2)
+  {
+    string variable = convertMatlabToString(prhs[i]);
+    if(variable == "ils_verbose")
+    {
+      string value = convertMatlabToString(prhs[i+1]);
+      if(value != "true" && value != "false")
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'ils_verbose\'. \'true\' or \'false\' expected.");
+      if(value == "true")
+        conf.sB("GPHIKClassifier", variable, true);
+      else
+        conf.sB("GPHIKClassifier", variable, false);
+    }
+
+    if(variable == "ils_max_iterations")
+    {
+      int value = convertMatlabToInt32(prhs[i+1]);
+      if(value < 1)
+        mexErrMsgIdAndTxt("mexnice:error","Expected parameter value larger than 0 for \'ils_max_iterations\'.");
+      conf.sI("GPHIKClassifier", variable, value);
+    }
+
+    if(variable == "ils_method")
+    {
+      string value = convertMatlabToString(prhs[i+1]);
+      if(value != "CG" && value != "CGL" && value != "SYMMLQ" && value != "MINRES")
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'ils_method\'. \'CG\', \'CGL\', \'SYMMLQ\' or \'MINRES\' expected.");
+        conf.sS("GPHIKClassifier", variable, value);
+    }
+
+    if(variable == "ils_min_delta")
+    {
+      double value = convertMatlabToDouble(prhs[i+1]);
+      if(value < 0.0)
+        mexErrMsgIdAndTxt("mexnice:error","Expected parameter value larger than 0 for \'ils_min_delta\'.");
+      conf.sD("GPHIKClassifier", variable, value);
+    }
+
+    if(variable == "ils_min_residual")
+    {
+      double value = convertMatlabToDouble(prhs[i+1]);
+      if(value < 0.0)
+        mexErrMsgIdAndTxt("mexnice:error","Expected parameter value larger than 0 for \'ils_min_residual\'.");
+      conf.sD("GPHIKClassifier", variable, value);
+    }
+
+
+    if(variable == "optimization_method")
+    {
+      string value = convertMatlabToString(prhs[i+1]);
+      if(value != "greedy" && value != "downhillsimplex" && value != "none")
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'optimization_method\'. \'greedy\', \'downhillsimplex\' or \'none\' expected.");
+        conf.sS("GPHIKClassifier", variable, value);
+    }
+
+    if(variable == "use_quantization")
+    {
+      string value = convertMatlabToString(prhs[i+1]);
+      if(value != "true" && value != "false")
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'use_quantization\'. \'true\' or \'false\' expected.");
+      if(value == "true")
+        conf.sB("GPHIKClassifier", variable, true);
+      else
+        conf.sB("GPHIKClassifier", variable, false);
+    }
+
+    if(variable == "num_bins")
+    {
+      int value = convertMatlabToInt32(prhs[i+1]);
+      if(value < 1)
+        mexErrMsgIdAndTxt("mexnice:error","Expected parameter value larger than 0 for \'num_bins\'.");
+      conf.sI("GPHIKClassifier", variable, value);
+    }
+
+    if(variable == "transform")
+    {
+      string value = convertMatlabToString(prhs[i+1]);
+      if(value != "absexp" && value != "exp" && value != "MKL" && value != "WeightedDim")
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'transform\'. \'absexp\', \'exp\' , \'MKL\' or \'WeightedDim\' expected.");
+        conf.sS("GPHIKClassifier", variable, value);
+    }
+
+    if(variable == "verboseTime")
+    {
+      string value = convertMatlabToString(prhs[i+1]);
+      if(value != "true" && value != "false")
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'verboseTime\'. \'true\' or \'false\' expected.");
+      if(value == "true")
+        conf.sB("GPHIKClassifier", variable, true);
+      else
+        conf.sB("GPHIKClassifier", variable, false);
+    }
+
+    if(variable == "verbose")
+    {
+      string value = convertMatlabToString(prhs[i+1]);
+      if(value != "true" && value != "false")
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'verbose\'. \'true\' or \'false\' expected.");
+      if(value == "true")
+        conf.sB("GPHIKClassifier", variable, true);
+      else
+        conf.sB("GPHIKClassifier", variable, false);
+    }
+
+    if(variable == "noise")
+    {
+      double value = convertMatlabToDouble(prhs[i+1]);
+      if(value < 0.0)
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value larger than 0 for \'noise\'.");
+      conf.sD("GPHIKClassifier", variable, value);
+    }
+
+    if(variable == "learn_balanced")
+    {
+      string value = convertMatlabToString(prhs[i+1]);
+      if(value != "true" && value != "false")
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'learn_balanced\'. \'true\' or \'false\' expected.");
+      if(value == "true")
+        conf.sB("GPHIKClassifier", variable, true);
+      else
+        conf.sB("GPHIKClassifier", variable, false);
+    }
+
+    if(variable == "optimize_noise")
+    {
+      string value = convertMatlabToString(prhs[i+1]);
+      if(value != "true" && value != "false")
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'optimize_noise\'. \'true\' or \'false\' expected.");
+      if(value == "true")
+        conf.sB("GPHIKClassifier", variable, true);
+      else
+        conf.sB("GPHIKClassifier", variable, false);
+    }
+    
+    if(variable == "varianceApproximation")
+    {
+      string value = convertMatlabToString(prhs[i+1]);
+      if(value != "approximate_fine" && value != "approximate_rough" && value != "exact" && value != "none")
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'varianceApproximation\'. \'approximate_fine\', \'approximate_rough\', \'none\' or \'exact\' expected.");
+        conf.sS("GPHIKClassifier", variable, value);
+    }
+    
+    if(variable == "nrOfEigenvaluesToConsiderForVarApprox")
+    {
+      double value = convertMatlabToDouble(prhs[i+1]);
+      conf.sI("GPHIKClassifier", variable, (int) value);
+    }    
+    
+  }
+
+
+  return conf;
+}
+
+// MAIN MATLAB FUNCTION
+void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
+{    
+    // get the command string specifying what to do
+    if (nrhs < 1)
+        mexErrMsgTxt("No commands and options passed... Aborting!");        
+    
+    if( !mxIsChar( prhs[0] ) )
+        mexErrMsgTxt("First argument needs to be the command, ie.e, the class method to call... Aborting!");        
+    
+    std::string cmd = convertMatlabToString( prhs[0] );
+      
+        
+    // create object
+    if ( !strcmp("new", cmd.c_str() ) )
+    {
+        // check output variable
+        if (nlhs != 1)
+            mexErrMsgTxt("New: One output expected.");
+        
+        // read config settings
+        NICE::Config conf = parseParameters(prhs+1,nrhs-1);
+        
+        // create class instance
+        NICE::GPHIKClassifier * classifier = new NICE::GPHIKClassifier ( &conf );
+        
+         
+        // handle to the C++ instance
+        plhs[0] = convertPtr2Mat<NICE::GPHIKClassifier>( classifier );
+        return;
+    }
+    
+    // in all other cases, there should be a second input,
+    // which the be the class instance handle
+    if (nrhs < 2)
+      mexErrMsgTxt("Second input should be a class instance handle.");
+    
+    // delete object
+    if ( !strcmp("delete", cmd.c_str() ) )
+    {
+        // Destroy the C++ object
+        destroyObject<NICE::GPHIKClassifier>(prhs[1]);
+        return;
+    }
+    
+    // get the class instance pointer from the second input
+    // every following function needs the classifier object
+    NICE::GPHIKClassifier * classifier = convertMat2Ptr<NICE::GPHIKClassifier>(prhs[1]);
+    
+    
+    ////////////////////////////////////////
+    //  Check which class method to call  //
+    ////////////////////////////////////////
+    
+    
+    // standard train - assumes initialized object
+    if (!strcmp("train", cmd.c_str() ))
+    {
+        // Check parameters
+        if (nlhs < 0 || nrhs < 4)
+        {
+            mexErrMsgTxt("Train: Unexpected arguments.");
+        }
+        
+        //------------- read the data --------------
+          
+        std::vector< NICE::SparseVector *> examplesTrain;
+        NICE::Vector yMultiTrain;  
+
+        if ( mxIsSparse( prhs[2] ) )
+        {
+            examplesTrain = convertSparseMatrixToNice( prhs[2] );
+        }
+        else
+        {
+            NICE::Matrix dataTrain;
+            dataTrain = convertMatrixToNice(prhs[2]);
+            
+            //----------------- convert data to sparse data structures ---------
+            examplesTrain.resize( dataTrain.rows() );
+
+                    
+            std::vector< NICE::SparseVector *>::iterator exTrainIt = examplesTrain.begin();
+            for (int i = 0; i < (int)dataTrain.rows(); i++, exTrainIt++)
+            {
+                *exTrainIt =  new NICE::SparseVector( dataTrain.getRow(i) );
+            }            
+        }
+          
+          yMultiTrain = convertVectorToNice(prhs[3]);
+          
+//           std::cerr << " DATA AFTER CONVERSION: \n" << std::endl;
+//           int lineIdx(0);
+//           for ( std::vector< NICE::SparseVector *>::const_iterator exTrainIt = examplesTrain.begin();
+//                 exTrainIt != examplesTrain.end(); exTrainIt++, lineIdx++)
+//           {
+//               std::cerr << "\n lineIdx: " << lineIdx << std::endl;
+//               (*exTrainIt)->store( std::cerr );
+//               
+//           }
+
+          // test assumption
+          {
+            if( yMultiTrain.Min() < 0)
+              mexErrMsgIdAndTxt("mexnice:error","Class labels smaller 0 are not allowed");
+          }
+
+
+          //----------------- train our classifier -------------
+          classifier->train ( examplesTrain , yMultiTrain );
+
+          //----------------- clean up -------------
+          for(int i=0;i<examplesTrain.size();i++)
+              delete examplesTrain[i];
+        
+        return;
+    }
+    
+    
+    // Classify    
+    if ( !strcmp("classify", cmd.c_str() ) )
+    {
+        // Check parameters
+        if ( (nlhs < 0) || (nrhs < 2) )
+        {
+            mexErrMsgTxt("Test: Unexpected arguments.");
+        }
+        
+        //------------- read the data --------------
+
+        int result;
+        NICE::SparseVector scores;
+        double uncertainty;        
+
+        if ( mxIsSparse( prhs[2] ) )
+        {
+            NICE::SparseVector * example;
+            example = new NICE::SparseVector ( convertSparseVectorToNice( prhs[2] ) );
+            classifier->classify ( example,  result, scores, uncertainty );
+            
+            //----------------- clean up -------------
+            delete example;
+        }
+        else
+        {
+            NICE::Vector * example;
+            example = new NICE::Vector ( convertVectorToNice(prhs[2]) ); 
+            classifier->classify ( example,  result, scores, uncertainty );
+            
+            //----------------- clean up -------------
+            delete example;            
+        }
+          
+          
+
+          // output
+          plhs[0] = mxCreateDoubleScalar( result ); 
+          
+          
+          if(nlhs >= 2)
+          {
+            plhs[1] = convertSparseVectorFromNice( scores, true  /*b_adaptIndex*/);
+          }
+          if(nlhs >= 3)
+          {
+            plhs[2] = mxCreateDoubleScalar( uncertainty );          
+          }
+          return;
+    }
+    
+    // Classify    
+    if ( !strcmp("uncertainty", cmd.c_str() ) )
+    {
+        // Check parameters
+        if ( (nlhs < 0) || (nrhs < 2) )
+        {
+            mexErrMsgTxt("Test: Unexpected arguments.");
+        }
+        
+        double uncertainty;        
+        
+        //------------- read the data --------------
+
+        if ( mxIsSparse( prhs[2] ) )
+        {
+            NICE::SparseVector * example;
+            example = new NICE::SparseVector ( convertSparseVectorToNice( prhs[2] ) );
+            classifier->predictUncertainty( example, uncertainty );
+            
+            //----------------- clean up -------------
+            delete example;            
+        }
+        else
+        {
+            NICE::Vector * example;
+            example = new NICE::Vector ( convertVectorToNice(prhs[2]) ); 
+            classifier->predictUncertainty( example, uncertainty );
+            
+            //----------------- clean up -------------
+            delete example;            
+        }
+        
+       
+
+          // output
+          plhs[0] = mxCreateDoubleScalar( uncertainty );                    
+          return;
+    }    
+    
+    
+    // Test    
+    if ( !strcmp("test", cmd.c_str() ) )
+    {        
+        // Check parameters
+        if (nlhs < 0 || nrhs < 4)
+            mexErrMsgTxt("Test: Unexpected arguments.");
+        //------------- read the data --------------
+        
+        
+        bool dataIsSparse ( mxIsSparse( prhs[2] ) );
+        
+        std::vector< NICE::SparseVector *> dataTest_sparse;
+        NICE::Matrix dataTest_dense;
+
+        if ( dataIsSparse )
+        {
+            dataTest_sparse = convertSparseMatrixToNice( prhs[2] );
+        }
+        else
+        {    
+            dataTest_dense = convertMatrixToNice(prhs[2]);          
+        }        
+
+          NICE::Vector yMultiTest;
+          yMultiTest = convertVectorToNice(prhs[3]);
+
+          
+          // ------------------------------------------
+          // ------------- PREPARATION --------------
+          // ------------------------------------------   
+          
+          // determine classes known during training and corresponding mapping
+          // thereby allow for non-continous class labels
+          std::set<int> classesKnownTraining = classifier->getKnownClassNumbers();
+          
+          int noClassesKnownTraining ( classesKnownTraining.size() );
+          std::map<int,int> mapClNoToIdxTrain;
+          std::set<int>::const_iterator clTrIt = classesKnownTraining.begin();
+          for ( int i=0; i < noClassesKnownTraining; i++, clTrIt++ )
+              mapClNoToIdxTrain.insert ( std::pair<int,int> ( *clTrIt, i )  );
+          
+          // determine classes known during testing and corresponding mapping
+          // thereby allow for non-continous class labels
+          std::set<int> classesKnownTest;
+          classesKnownTest.clear();
+          
+  
+          // determine which classes we have in our label vector
+          // -> MATLAB: myClasses = unique(y);
+          for ( NICE::Vector::const_iterator it = yMultiTest.begin(); it != yMultiTest.end(); it++ )
+          {
+            if ( classesKnownTest.find ( *it ) == classesKnownTest.end() )
+            {
+              classesKnownTest.insert ( *it );
+            }
+          }          
+          
+          int noClassesKnownTest ( classesKnownTest.size() );  
+          std::map<int,int> mapClNoToIdxTest;
+          std::set<int>::const_iterator clTestIt = classesKnownTest.begin();
+          for ( int i=0; i < noClassesKnownTest; i++, clTestIt++ )
+              mapClNoToIdxTest.insert ( std::pair<int,int> ( *clTestIt, i )  );          
+          
+
+
+          int i_numTestSamples;
+          
+          if ( dataIsSparse ) 
+              i_numTestSamples = dataTest_sparse.size();
+          else
+              i_numTestSamples = (int) dataTest_dense.rows();
+          
+          NICE::Matrix confusionMatrix( noClassesKnownTraining, noClassesKnownTest, 0.0);
+          NICE::Matrix scores( i_numTestSamples, noClassesKnownTraining, 0.0);
+          
+          
+
+          // ------------------------------------------
+          // ------------- CLASSIFICATION --------------
+          // ------------------------------------------          
+          
+          NICE::Timer t;
+          double testTime (0.0);
+          
+
+
+          for (int i = 0; i < i_numTestSamples; i++)
+          {
+             //----------------- convert data to sparse data structures ---------
+            
+
+             int result;
+             NICE::SparseVector exampleScoresSparse;
+
+             if ( dataIsSparse )
+             {                
+                // and classify
+                t.start();
+                classifier->classify( dataTest_sparse[ i ], result, exampleScoresSparse );
+                t.stop();
+                testTime += t.getLast();
+             }
+             else
+             {
+                 NICE::Vector example ( dataTest_dense.getRow(i) );
+                // and classify
+                t.start();
+                classifier->classify( &example, result, exampleScoresSparse );
+                t.stop();
+                testTime += t.getLast();                
+             }
+
+             confusionMatrix(  mapClNoToIdxTrain.find(result)->second, mapClNoToIdxTest.find(yMultiTest[i])->second ) += 1.0;
+             int scoreCnt ( 0 );
+             for ( NICE::SparseVector::const_iterator scoreIt = exampleScoresSparse.begin(); scoreIt != exampleScoresSparse.end(); scoreIt++, scoreCnt++ )
+                scores(i,scoreCnt) = scoreIt->second;
+                
+          }
+          
+          std::cerr << "Time for testing: " << testTime << std::endl;          
+          
+          // clean up
+          if ( dataIsSparse )
+          {
+              for ( std::vector<NICE::SparseVector *>::iterator it = dataTest_sparse.begin(); it != dataTest_sparse.end(); it++) 
+                  delete *it;
+          }
+          
+
+
+          confusionMatrix.normalizeColumnsL1();
+          //std::cerr << confusionMatrix << std::endl;
+
+          double recRate = confusionMatrix.trace()/confusionMatrix.rows();
+          //std::cerr << "average recognition rate: " << recRate << std::endl;
+
+          
+          plhs[0] = mxCreateDoubleScalar( recRate );
+
+          if(nlhs >= 2)
+            plhs[1] = convertMatrixFromNice(confusionMatrix);
+          if(nlhs >= 3)
+            plhs[2] = convertMatrixFromNice(scores);          
+          
+          
+        return;
+    }
+    
+    // store the classifier    
+    if ( !strcmp("store", cmd.c_str() ) || !strcmp("save", cmd.c_str() ) )
+    {
+        // Check parameters
+        if ( nrhs < 3 )
+            mexErrMsgTxt("store: no destination given.");        
+               
+        std::string s_destination = convertMatlabToString( prhs[2] );
+          
+        std::filebuf fb;
+        fb.open ( s_destination.c_str(), ios::out );
+        std::ostream os(&fb);
+        //
+        classifier->store( os );
+        //   
+        fb.close();        
+            
+        return;
+    }
+    
+    // load classifier from external file    
+    if ( !strcmp("restore", cmd.c_str() ) || !strcmp("load", cmd.c_str() ) )
+    {
+        // Check parameters
+        if ( nrhs < 3 )
+            mexErrMsgTxt("restore: no destination given.");        
+               
+        std::string s_destination = convertMatlabToString( prhs[2] );
+          
+        std::filebuf fbIn;
+        fbIn.open ( s_destination.c_str(), ios::in );
+        std::istream is (&fbIn);
+        //
+        classifier->restore( is );
+        //   
+        fbIn.close();        
+            
+        return;
+    }    
+    
+    
+    // Got here, so command not recognized
+    
+    std::string errorMsg (cmd.c_str() );
+    errorMsg += " -- command not recognized.";
+    mexErrMsgTxt( errorMsg.c_str() );
+
+}

matlab/Makefile

@@ -0,0 +1,5 @@
+NICEFLAGS1=$(shell pkg-config libgp-hik-core --cflags --libs)
+NICEFLAGS=$(subst -fopenmp,,$(NICEFLAGS1))
+
+default:
+	/home/matlab/7.14/bin/mex ${NICEFLAGS} -largeArrayDims GPHIK.cpp 

matlab/classHandleMtoC.h

@@ -0,0 +1,142 @@
+/** 
+* @file classHandleMtoC.h
+* @brief Generic class to pass C++ objects to matlab (Interface and inline implementations)
+* @author Alexander Freytag
+* @date 19-12-2013 (dd-mm-yyyy)
+
+*/
+#ifndef _NICE_CLASSHANDLEMTOCINCLUDE
+#define _NICE_CLASSHANDLEMTOCINCLUDE
+
+#include "mex.h"
+#include <stdint.h>
+#include <iostream>
+#include <string>
+#include <cstring>
+#include <typeinfo>
+
+#define CLASS_HANDLE_SIGNATURE 0xFF00F0A3
+
+  /** 
+  * @class FMKGPHyperparameterOptimization
+  * @brief Generic class to pass C++ objects to matlab
+  * @author Alexander Freytag
+  */
+template<class objectClass> class ClassHandle
+{
+  private:
+      //!
+      uint32_t i_mySignature;
+      
+      //! typeid.name of object class we refere to
+      std::string s_myName;
+      
+      //! the actual pointer to our C++ object
+      objectClass* p_myPtr;  
+  
+  public:
+    
+    /**
+    * @brief standard constructor
+    *
+    * @param ptr pointer to the c++ object
+    */    
+      ClassHandle ( objectClass* p_ptr ) : p_myPtr(p_ptr), s_myName( typeid(objectClass).name() )
+      {
+        i_mySignature = CLASS_HANDLE_SIGNATURE;
+      }
+      
+    /**
+    * @brief standard destructor
+    */        
+      ~ClassHandle()
+      {
+          // reset internal variables
+          i_mySignature = 0;
+          
+          // clearn up data
+          delete p_myPtr;
+      }
+      
+    /**
+    * @brief check whether the class handle was initialized properly, i.e., we point to an actual object
+    */       
+      bool isValid()
+      { 
+        return ( (i_mySignature == CLASS_HANDLE_SIGNATURE) && !strcmp( s_myName.c_str(), typeid(objectClass).name() )   );
+      }
+      
+    /**
+    * @brief get the pointer to the actual object
+    */       
+      objectClass * getPtrToObject()
+      { 
+        return p_myPtr;
+      }
+
+
+};
+
+
+////////////////////////////////////////////
+//           conversion methods           //
+////////////////////////////////////////////
+
+/**
+* @brief convert handle to C++ object into matlab usable data
+*/ 
+template<class objectClass> inline mxArray *convertPtr2Mat(objectClass *ptr)
+{
+    // prevent user from clearing the mex file! Otherwise, storage leaks might be caused
+    mexLock();
+    
+    // allocate memory
+    mxArray *out = mxCreateNumericMatrix(1, 1, mxUINT64_CLASS, mxREAL);
+    
+    // convert handle do matlab usable data
+    *((uint64_t *)mxGetData(out)) = reinterpret_cast<uint64_t>(new ClassHandle<objectClass>(ptr));
+    
+    return out;
+}
+
+/**
+* @brief convert matlab usable data referring to an object into handle to C++ object
+*/ 
+template<class objectClass> inline ClassHandle<objectClass> *convertMat2HandlePtr(const mxArray *in)
+{
+    // check that the given pointer actually points to a real object
+    if ( ( mxGetNumberOfElements(in) != 1 )     ||
+         ( mxGetClassID(in) != mxUINT64_CLASS ) ||
+           mxIsComplex(in)
+       )
+        mexErrMsgTxt("Input must be a real uint64 scalar.");
+        
+    ClassHandle<objectClass> *ptr = reinterpret_cast<ClassHandle<objectClass> *>(*((uint64_t *)mxGetData(in)));
+    
+    if (!ptr->isValid())
+        mexErrMsgTxt("Handle not valid.");
+    
+    return ptr;
+}
+
+/**
+* @brief convert matlab usable data referring to an object into direct pointer to the underlying C++ object
+*/ 
+template<class objectClass> inline objectClass *convertMat2Ptr(const mxArray *in)
+{
+    return convertMat2HandlePtr<objectClass>(in)->getPtrToObject();
+}
+
+/**
+* @brief convert matlab usable data referring to an object into direct pointer to the underlying C++ object
+*/ 
+template<class objectClass> inline void destroyObject(const mxArray *in)
+{
+    // clean up
+    delete convertMat2HandlePtr<objectClass>(in);
+    
+    // storage is freed, so users can savely clear the mex file again at any time...
+    mexUnlock();
+}
+
+#endif // _NICE_CLASSHANDLEMTOCINCLUDE

progs/toyExample.cpp

@@ -23,6 +23,7 @@ int main (int argc, char* argv[])
   
   Config conf ( argc, argv );
   std::string trainData = conf.gS( "main", "trainData", "progs/toyExampleSmallScaleTrain.data" );
+  bool b_debug = conf.gB( "main", "debug", false );
 
   
   //------------- read the training data --------------
@@ -31,19 +32,38 @@ int main (int argc, char* argv[])
   NICE::Vector yBinTrain;
   NICE::Vector yMultiTrain;  
 
-  std::ifstream ifsTrain ( trainData.c_str() , ios::in );
-
-  if (ifsTrain.good() )
-  {
-    ifsTrain >> dataTrain;
-    ifsTrain >> yBinTrain;
-    ifsTrain >> yMultiTrain;
-    ifsTrain.close();  
+  if ( b_debug )
+  { 
+    dataTrain.resize(6,3);
+    dataTrain.set(0);
+    dataTrain(0,0) = 0.2; dataTrain(0,1) = 0.3; dataTrain(0,2) = 0.5;
+    dataTrain(1,0) = 0.3; dataTrain(1,1) = 0.2; dataTrain(1,2) = 0.5;    
+    dataTrain(2,0) = 0.9; dataTrain(2,1) = 0.0; dataTrain(2,2) = 0.1;
+    dataTrain(3,0) = 0.8; dataTrain(3,1) = 0.1; dataTrain(3,2) = 0.1;    
+    dataTrain(4,0) = 0.1; dataTrain(4,1) = 0.1; dataTrain(4,2) = 0.8;
+    dataTrain(5,0) = 0.1; dataTrain(5,1) = 0.0; dataTrain(5,2) = 0.9;    
+    
+    yMultiTrain.resize(6);
+    yMultiTrain[0] = 1; yMultiTrain[1] = 1;
+    yMultiTrain[2] = 2; yMultiTrain[3] = 2;
+    yMultiTrain[2] = 3; yMultiTrain[3] = 3;
   }
   else 
   {
-    std::cerr << "Unable to read training data, aborting." << std::endl;
-    return -1;
+    std::ifstream ifsTrain ( trainData.c_str() , ios::in );
+
+    if (ifsTrain.good() )
+    {
+      ifsTrain >> dataTrain;
+      ifsTrain >> yBinTrain;
+      ifsTrain >> yMultiTrain;
+      ifsTrain.close();  
+    }
+    else 
+    {
+      std::cerr << "Unable to read training data, aborting." << std::endl;
+      return -1;
+    }
   }
   
   //----------------- convert data to sparse data structures ---------
@@ -56,6 +76,8 @@ int main (int argc, char* argv[])
     *exTrainIt =  new NICE::SparseVector( dataTrain.getRow(i) );
   }
   
+  std::cerr << "Number of training examples: " << examplesTrain.size() << std::endl;
+  
   //----------------- train our classifier -------------
   conf.sB("GPHIKClassifier", "verbose", false);
   GPHIKClassifier * classifier  = new GPHIKClassifier ( &conf );  
@@ -68,23 +90,36 @@ int main (int argc, char* argv[])
   
   //------------- read the test data --------------
   
+  
   NICE::Matrix dataTest;
   NICE::Vector yBinTest;
-  NICE::Vector yMultiTest;  
-  
-  std::string testData = conf.gS( "main", "testData", "progs/toyExampleTest.data" );  
-  std::ifstream ifsTest ( testData.c_str(), ios::in );
-  if (ifsTest.good() )
-  {
-    ifsTest >> dataTest;
-    ifsTest >> yBinTest;
-    ifsTest >> yMultiTest;
-    ifsTest.close();  
+  NICE::Vector yMultiTest; 
+    
+  if ( b_debug )
+  { 
+    dataTest.resize(1,3);
+    dataTest.set(0);
+    dataTest(0,0) = 0.3; dataTest(0,1) = 0.4; dataTest(0,2) = 0.3;
+    
+    yMultiTrain.resize(1);
+    yMultiTrain[0] = 1;
   }
   else 
-  {
-    std::cerr << "Unable to read test data, aborting." << std::endl;
-    return -1;
+  {  
+    std::string testData = conf.gS( "main", "testData", "progs/toyExampleTest.data" );  
+    std::ifstream ifsTest ( testData.c_str(), ios::in );
+    if (ifsTest.good() )
+    {
+      ifsTest >> dataTest;
+      ifsTest >> yBinTest;
+      ifsTest >> yMultiTest;
+      ifsTest.close();  
+    }
+    else 
+    {
+      std::cerr << "Unable to read test data, aborting." << std::endl;
+      return -1;
+    }
   }
   
   //TODO adapt this to the actual number of classes
@@ -93,7 +128,16 @@ int main (int argc, char* argv[])
   NICE::Timer t;
   double testTime (0.0);
   
-  for (int i = 0; i < (int)dataTest.rows(); i++)
+  double uncertainty;
+  
+  int i_loopEnd  ( (int)dataTest.rows() );
+  
+  if ( b_debug )
+  {
+    i_loopEnd = 1;
+  }
+  
+  for (int i = 0; i < i_loopEnd ; i++)
   {
     //----------------- convert data to sparse data structures ---------
     NICE::SparseVector * example =  new NICE::SparseVector( dataTest.getRow(i) );
@@ -107,15 +151,29 @@ int main (int argc, char* argv[])
     t.stop();
     testTime += t.getLast();
     
-    confusionMatrix(result, yMultiTest[i]) += 1.0;
+    std::cerr << " scores.size(): " << scores.size() << std::endl;
+    scores.store(std::cerr);
+    
+    if ( b_debug )
+    {    
+      classifier->predictUncertainty( example, uncertainty );
+      std::cerr << " uncertainty: " << uncertainty << std::endl;
+    }
+    else
+    {
+      confusionMatrix(result, yMultiTest[i]) += 1.0;
+    }
   }
   
-  std::cerr << "Time for testing: " << testTime << std::endl;
-  
-  confusionMatrix.normalizeColumnsL1();
-  std::cerr << confusionMatrix << std::endl;
+  if ( !b_debug )
+  {
+    std::cerr << "Time for testing: " << testTime << std::endl;
+    
+    confusionMatrix.normalizeColumnsL1();
+    std::cerr << confusionMatrix << std::endl;
 
-  std::cerr << "average recognition rate: " << confusionMatrix.trace()/confusionMatrix.rows() << std::endl;
+    std::cerr << "average recognition rate: " << confusionMatrix.trace()/confusionMatrix.rows() << std::endl;
+  }
   
   
   return 0;

tests/TestFastHIK.cpp

@@ -15,28 +15,6 @@
 
 #include "TestFastHIK.h"
 
-
-const bool verbose = false;
-const bool verboseStartEnd = true;
-const bool solveLinWithoutRand = false;
-const uint n = 20;//1500;//1500;//10;
-const uint d = 5;//200;//2;
-const uint numBins = 11;//1001;//1001;
-const uint solveLinMaxIterations = 1000;
-const double sparse_prob = 0.6;
-const bool smallTest = false;
-
-using namespace NICE;
-using namespace std;
-
-CPPUNIT_TEST_SUITE_REGISTRATION( TestFastHIK );
-
-void TestFastHIK::setUp() {
-}
-
-void TestFastHIK::tearDown() {
-}
-
 bool compareVVector(const NICE::VVector & A, const NICE::VVector & B, const double & tolerance = 10e-8)
 {
   bool result(true);
@@ -91,6 +69,27 @@ bool compareLUTs(const double* LUT1, const double* LUT2, const int & size, const
   return result;
 }
 
+const bool verbose = false;
+const bool verboseStartEnd = true;
+const bool solveLinWithoutRand = false;
+const uint n = 20;//1500;//1500;//10;
+const uint d = 5;//200;//2;
+const uint numBins = 11;//1001;//1001;
+const uint solveLinMaxIterations = 1000;
+const double sparse_prob = 0.6;
+const bool smallTest = false;
+
+using namespace NICE;
+using namespace std;
+
+CPPUNIT_TEST_SUITE_REGISTRATION( TestFastHIK );
+
+void TestFastHIK::setUp() {
+}
+
+void TestFastHIK::tearDown() {
+}
+
 void TestFastHIK::testKernelMultiplication() 
 {
   if (verboseStartEnd)
@@ -789,454 +788,4 @@ void TestFastHIK::testKernelVector()
   
 }
 
-void TestFastHIK::testAddExample()
-{
-  if (verboseStartEnd)
-    std::cerr << "================== TestFastHIK::testAddExample ===================== " << std::endl;  
-  
-  std::vector< std::vector<double> > dataMatrix;
-  int dim = 3;
-  int number = 5;
-  
-  if (!smallTest)
-  {
-    dim = d;
-    number = n;
-  }
-  
-  if (smallTest)
-  {
-    dataMatrix.resize(3);
-    //we explicitely give some values which can easily be verified
-    dataMatrix[0].push_back(0.2);dataMatrix[0].push_back(0.1);dataMatrix[0].push_back(0.0);dataMatrix[0].push_back(0.0);dataMatrix[0].push_back(0.4); 
-    dataMatrix[1].push_back(0.3);dataMatrix[1].push_back(0.6);dataMatrix[1].push_back(1.0);dataMatrix[1].push_back(0.4);dataMatrix[1].push_back(0.3);
-    dataMatrix[2].push_back(0.5);dataMatrix[2].push_back(0.3);dataMatrix[2].push_back(0.0);dataMatrix[2].push_back(0.6);dataMatrix[2].push_back(0.3);
-  }
-  else
-  {
-    // randomly generate features
-    generateRandomFeatures ( dim, number, dataMatrix );
-
-    // and make them sparse
-    int nrZeros(0);
-    for ( int i = 0 ; i < dim; i++ )
-    {
-      for ( int k = 0; k < number; k++ )
-        if ( drand48() < sparse_prob ) 
-        {
-          dataMatrix[i][k] = 0.0;
-          nrZeros++;
-        }
-    }    
-  }
-  
-  if ( verbose ) {
-    std::cerr << "data matrix: " << std::endl;
-    printMatrix ( dataMatrix );
-    std::cerr << endl;
-  }
-  
-  double noise = 1.0;
-  //check the features stored in the fmk
-  FastMinKernel fmk ( dataMatrix, noise );  
-  NICE::Vector alpha;
-  
-  ParameterizedFunction *pf = new PFAbsExp( 1.2 ); //1.0 is okay
-  fmk.applyFunctionToFeatureMatrix( pf );
-//   pf->applyFunctionToFeatureMatrix ( fmk.featureMatrix() );  
-  
-  std::cerr << "generate alpha" << std::endl;
-  
-  if (smallTest)
-  {
-    //we explicitely give some values which can easily be verified
-    alpha = Vector(5,1.0);
-    alpha[0] = 0.1;alpha[1] = 0.2;alpha[2] = 0.4;alpha[3] = 0.8;alpha[4] = 1.6;
-  }
-  else
-  {  // randomly generate features
-     alpha = Vector::UniformRandom( number, 0.0, 1.0, 0 );
-  }
-  
-  
-  std::cerr << "generate xStar" << std::endl;
-  std::vector<double> xStar;
-  if (smallTest)
-  {
-    // we check the following cases: largest elem in dim, smallest elem in dim, zero element
-    // remember to adapt the feature in some lines apart as well    
-    xStar.push_back(0.9);xStar.push_back(0.0);xStar.push_back(0.1);
-  }
-  else
-  {
-    // again: random sampling
-    for ( int i = 0 ; i < dim; i++ )
-    {
-      if ( drand48() < sparse_prob ) 
-        xStar.push_back(0.0);
-      else
-        xStar.push_back(drand48());
-    }
-  }
-  NICE::Vector xStarVec (xStar);
-  NICE::SparseVector xStarSV (xStarVec);
-  
-  // check the alpha-preparations
-  NICE::VVector A;
-  NICE::VVector B;
-  fmk.hik_prepare_alpha_multiplications( alpha, A, B );
-  
-  //check the quantization and LUT construction
-  Quantization q ( numBins );  
-  //direct
-//   double * LUT = fmk.hikPrepareLookupTable(alpha, q);
-  //indirect
-  double * LUT = fmk.hik_prepare_alpha_multiplications_fast( A, B, q, pf );
-  
-  //check for kernel vector norm approximation
-  NICE::VVector AForKVN;
-  fmk.hikPrepareKVNApproximation(AForKVN);
-  
-  //check the LUTs for fast kernel vector norm approximation
-  //direct
-  double* LUT_kernelVectorNormDirect = fmk.hikPrepareLookupTableForKVNApproximation(q, pf );
-  //indirect
-  double* LUT_kernelVectorNorm = fmk.hikPrepareKVNApproximationFast( AForKVN, q, pf );
-  
-  bool LUTKVN_equal( compareLUTs( LUT_kernelVectorNorm, LUT_kernelVectorNormDirect, q.size()*dim ) );
-  
-  if (verbose)
-  {
-    if (LUTKVN_equal == false)
-    {
-      std::cerr << "LUTKVN is not equal :( " << std::endl;
-        std::cerr << "LUT_kernelVectorNorm: " << std::endl;
-        for ( uint i = 0; i < q.size()*dim; i++ )
-        {
-          if ( (i % q.size()) == 0)
-            std::cerr << std::endl;
-          std::cerr << LUT_kernelVectorNorm[i] << " ";
-        }
-        std::cerr << "LUT_kernelVectorNormDirect: "<< std::endl;
-        for ( uint i = 0; i < q.size()*dim; i++ )
-        {
-          if ( (i % q.size()) == 0)
-            std::cerr << std::endl;
-          std::cerr << LUT_kernelVectorNormDirect[i] << " ";
-        }      
-    }
-  }
-  CPPUNIT_ASSERT( LUTKVN_equal == true );
-  
-  if (verbose)
-    std::cerr << "start the incremental learning part" << std::endl;
-
-  // ------  Incremental Learning -----
-  
-  double newAlpha;
-  if (smallTest) 
-    newAlpha = 3.2;
-  else
-    newAlpha = drand48();
-  alpha.append(newAlpha);
-   
-  // add an example
-  if (verbose)
-    std::cerr << "addExample" << std::endl;  
-  fmk.addExample( xStarSV, pf );  
-  
-  // update the alpha preparation
-  if (verbose)  
-    std::cerr << "update Alpha Preparation" << std::endl;
-  fmk.updatePreparationForAlphaMultiplications( xStarSV, newAlpha, A, B, pf );
-  
-  // update the LUT for fast multiplications
-  if (verbose)  
-    std::cerr << "update LUT" << std::endl;
-  fmk.updateLookupTableForAlphaMultiplications( xStarSV, newAlpha, LUT, q, pf );
-  
-  //update VVector for Kernel vector norm
-  if (verbose)  
-    std::cerr << "update VVector for Kernel vector norm" << std::endl;
-  fmk.updatePreparationForKVNApproximation( xStarSV, AForKVN, pf );
-  
-  // update LUT for kernel vector norm
-  if (verbose)  
-    std::cerr << "update LUT for kernel vector norm" << std::endl;
-  fmk.updateLookupTableForKVNApproximation( xStarSV, LUT_kernelVectorNorm, q, pf );
-  
-  //and batch retraining  
-  if (verbose)  
-    std::cerr << "perform batch retraining " << std::endl;  
-  for ( int i = 0 ; i < dim; i++ )
-    dataMatrix[i].push_back(xStar[i]);
-  
-  FastMinKernel fmk2 ( dataMatrix, noise );
-  fmk2.applyFunctionToFeatureMatrix( pf );
-  
-  NICE::VVector A2;
-  NICE::VVector B2;
-  if (verbose)  
-    std::cerr << "prepare alpha multiplications" << std::endl;
-  fmk2.hik_prepare_alpha_multiplications( alpha, A2, B2 );
- 
-  // compare the content of the data matrix
-  if (verbose)  
-    std::cerr << "do the comparison of the resulting feature matrices" << std::endl;
-  if (verbose)
-  {
-    std::cerr << "fmk.featureMatrix().print()" << std::endl;
-    fmk.featureMatrix().print(std::cerr);
-  
-    std::cerr << "fmk2.featureMatrix().print()" << std::endl;
-    fmk2.featureMatrix().print(std::cerr);
-  }  
-  
-  CPPUNIT_ASSERT(fmk.featureMatrix() == fmk2.featureMatrix());
-
-  //compare the preparation for alpha multiplications
-  if (verbose)  
-    std::cerr << "do the comparison of the resulting matrices A and B" << std::endl;
-  CPPUNIT_ASSERT(compareVVector(A, A2));  
-  CPPUNIT_ASSERT(compareVVector(B, B2));
-  
-  if (verbose)
-  {
-    std::cerr << "compare the preparation for alpha multiplications" << std::endl;
-    std::cerr << "A: " << std::endl;
-    A.store(std::cerr);
-    std::cerr << "A2: " << std::endl;
-    A2.store(std::cerr);
-    std::cerr << "B: " << std::endl;
-    B.store(std::cerr);
-    std::cerr << "B2: " << std::endl;
-    B2.store(std::cerr);
-  }  
-  
-  // compare the resulting LUTs
-  if (verbose)
-    std::cerr << "prepare LUT" << std::endl;
-  double * LUT2 = fmk2.hikPrepareLookupTable( alpha, q, pf );    
-  if (verbose)
-    std::cerr << "do the comparison of the resulting LUTs" << std::endl;  
-  bool LUTequal( compareLUTs( LUT, LUT2, q.size()*dim) );
-  
-  if (verbose)
-  {
-    if ( LUTequal )
-      std::cerr << "LUTs are equal :) " << std::endl;
-    else
-    {
-      std::cerr << "LUTs are not equal :( " << std::endl;
-      std::cerr << "new feature vector: " << xStarVec << std::endl;
-      
-      std::cerr << "newAlpha: " << newAlpha <<  " alpha " << alpha << std::endl;
-      std::cerr << "LUT: " << std::endl;
-      for ( uint i = 0; i < q.size()*dim; i++ )
-      {
-        if ( (i % q.size()) == 0)
-          std::cerr << std::endl;
-        std::cerr << LUT[i] << " ";
-      }
-      std::cerr << "LUT2: "<< std::endl;
-      for ( uint i = 0; i < q.size()*dim; i++ )
-      {
-        if ( (i % q.size()) == 0)
-          std::cerr << std::endl;
-        std::cerr << LUT2[i] << " ";
-      }     
-    }
-  }
-  CPPUNIT_ASSERT( LUTequal );
-  
-  //check for kernel vector norm approximation
-  NICE::VVector A2ForKVN;
-  fmk2.hikPrepareKVNApproximation( A2ForKVN );
-  bool KVN_equal ( compareVVector(AForKVN, A2ForKVN) );
- 
-  if (verbose)
-  {
-    if ( KVN_equal )
-      std::cerr << "VVectors for kernel vector norm are equal :) " << std::endl;
-    else
-    {
-      std::cerr << "VVectors for vector norm are not equal :( " << std::endl;
-      std::cerr << "new feature vector: " << xStarVec << std::endl;
-      
-      std::cerr << "AForKVN: " << std::endl;
-      AForKVN.store(std::cerr);
-      
-      std::cerr << "A2ForKVN: "<< std::endl;
-      A2ForKVN.store(std::cerr);
-    }
-  }  
-  
-  CPPUNIT_ASSERT( KVN_equal );  
-  
-  //check for kernel vector norm approximation with LUTs
-  if (verbose)
-    std::cerr << "prepare LUT for kernel vector norm" << std::endl;
-  double* LUT2_kernelVectorNorm = fmk2.hikPrepareLookupTableForKVNApproximation( q, pf );  
-  if (verbose)
-    std::cerr << "do the comparison of the resulting LUTs for kernel vector norm computation" << std::endl;
-  bool LUT_KVN_equal( compareLUTs ( LUT_kernelVectorNorm, LUT2_kernelVectorNorm, q.size()*dim ) );
-  
-  if (verbose)
-  {
-    if ( LUT_KVN_equal )
-      std::cerr << "LUTs for kernel vector norm are equal :) " << std::endl;
-    else
-    {
-      std::cerr << "LUTs kernel vector norm are not equal :( " << std::endl;
-      std::cerr << "new feature vector: " << xStarVec << std::endl;
-      
-      std::cerr << "LUT_kernelVectorNorm: " << std::endl;
-      for ( int i = 0; i < q.size()*dim; i++ )
-      {
-        if ( (i % q.size()) == 0)
-          std::cerr << std::endl;
-        std::cerr << LUT_kernelVectorNorm[i] << " ";
-      }
-      std::cerr << std::endl << "LUT2_kernelVectorNorm: "<< std::endl;
-      for ( uint i = 0; i < q.size()*dim; i++ )
-      {
-        if ( (i % q.size()) == 0)
-          std::cerr << std::endl;
-        std::cerr << LUT2_kernelVectorNorm[i] << " ";
-      }     
-    }
-  }  
-  
-  CPPUNIT_ASSERT( LUT_KVN_equal );
-  
-  delete [] LUT;
-  delete [] LUT2;
-  
-  delete [] LUT_kernelVectorNorm;
-  delete [] LUT2_kernelVectorNorm;
-  
-  if (verboseStartEnd)
-    std::cerr << "================== TestFastHIK::testAddExample done ===================== " << std::endl;  
-}
-
-void TestFastHIK::testAddMultipleExamples()
-{
-  if (verboseStartEnd)
-    std::cerr << "================== TestFastHIK::testAddMultipleExamples ===================== " << std::endl;  
-  
-  std::vector< std::vector<double> > dataMatrix;
-  int dim = d;
-  int number = n;
-
-  // randomly generate features
-  generateRandomFeatures ( dim, number, dataMatrix );
-
-  // and make them sparse
-  int nrZeros(0);
-  for ( int i = 0 ; i < dim; i++ )
-  {
-    for ( int k = 0; k < number; k++ )
-      if ( drand48() < sparse_prob ) 
-      {
-        dataMatrix[i][k] = 0.0;
-        nrZeros++;
-      }
-  }    
-  
-  if ( verbose ) {
-    std::cerr << "data matrix: " << std::endl;
-    printMatrix ( dataMatrix );
-    std::cerr << endl;
-  }
-  
-  double noise = 1.0;
-  //check the features stored in the fmk
-  FastMinKernel fmk ( dataMatrix, noise );  
-  NICE::Vector alpha;
-  
-  ParameterizedFunction *pf = new PFAbsExp( 1.0 ); //1.0 is okay
-  fmk.applyFunctionToFeatureMatrix( pf );
-  
-  std::cerr << "generate alpha" << std::endl;  
-  
-  // randomly generate features
-  alpha = Vector::UniformRandom( number, 0.0, 1.0, 0 );
-   
-/*  // check the alpha-preparations
-  NICE::VVector A;
-  NICE::VVector B;
-  fmk.hik_prepare_alpha_multiplications( alpha, A, B );*/  
-  
-  if (verbose)
-    std::cerr << "start the incremental learning part" << std::endl;
-
-  // ------  Incremental Learning -----
-    
-  std::cerr << "generate xStar" << std::endl;
-  std::vector<NICE::SparseVector > newExamples;
-  int nrOfNewExamples(5);
-  // again: random sampling
-  for (int i = 0; i < nrOfNewExamples; i++)
-  {
-    NICE::Vector xStar(dim);
-    for ( int j = 0 ; j < dim; j++ )
-    {
-      if ( drand48() < sparse_prob ) 
-      {
-        xStar[j] = 0.0;
-        dataMatrix[j].push_back(0.0);
-      }
-      else
-      {
-        double tmp(drand48());
-        xStar[j] = tmp;
-        dataMatrix[j].push_back(tmp);
-      }
-    }
-    
-    NICE::SparseVector xStarSV (xStar);
-    newExamples.push_back(xStarSV);
-  }    
-
-  // add an example
-  if (verbose)
-    std::cerr << "addExample" << std::endl;  
-  for (int i = 0; i < nrOfNewExamples; i++)
-  {
-    fmk.addExample( newExamples[i], pf );  
-  }
-  
-  int oldSize(alpha.size());
-  alpha.resize( oldSize + nrOfNewExamples);
-  for (int i = 0; i < nrOfNewExamples; i++)
-  {
-    alpha[oldSize + i] = drand48();
-  }
-   
-  
-  // update the alpha preparation
-  if (verbose)  
-    std::cerr << "update Alpha Preparation" << std::endl;
-  // check the alpha-preparations
-  NICE::VVector A;
-  NICE::VVector B;
-  fmk.hik_prepare_alpha_multiplications( alpha, A, B );   
-  
-  FastMinKernel fmk2 ( dataMatrix, noise );
-  fmk2.applyFunctionToFeatureMatrix( pf );  
-  
-  NICE::VVector A2;
-  NICE::VVector B2;
-  fmk2.hik_prepare_alpha_multiplications( alpha, A2, B2 );
-  
-  bool equalA = compareVVector( A, A2 );
-  bool equalB = compareVVector( B, B2 );
-  
-  CPPUNIT_ASSERT(equalA == true);
-  CPPUNIT_ASSERT(equalB == true);  
-  
-  if (verboseStartEnd)
-    std::cerr << "================== TestFastHIK::testAddMultipleExamples done ===================== " << std::endl;  
-}
-
 #endif

tests/TestFastHIK.h

@@ -19,8 +19,6 @@ class TestFastHIK : public CppUnit::TestFixture {
     CPPUNIT_TEST(testLUTUpdate);
     CPPUNIT_TEST(testLinSolve);
     CPPUNIT_TEST(testKernelVector);
-    CPPUNIT_TEST(testAddExample);
-    CPPUNIT_TEST(testAddMultipleExamples);
     
     CPPUNIT_TEST_SUITE_END();
   
@@ -40,8 +38,6 @@ class TestFastHIK : public CppUnit::TestFixture {
     void testLUTUpdate();
     void testLinSolve();
     void testKernelVector();
-    void testAddExample();
-    void testAddMultipleExamples();
 
 };