first but broken direct implementation of the partial matrix idea

FastMinKernel.cpp

@@ -227,6 +227,59 @@ void FastMinKernel::hik_prepare_alpha_multiplications(const NICE::Vector & _alph
 
 }
 
+void FastMinKernel::hik_prepare_alpha_multiplications_subset(const NICE::Vector & _alpha,
+                                                      NICE::VVector & _A,
+                                                      NICE::VVector & _B,
+                                                      const vector<int> & _subset) const
+{
+  _A.resize( this->ui_d );
+  _B.resize( this->ui_d );
+
+  //  we only need as many entries as we have nonZero entries in our features for the corresponding dimensions
+  for (uint i = 0; i < this->ui_d; i++)
+  {
+    uint numNonZero = this->X_sorted.getNumberOfNonZeroElementsPerDimension(i);
+    _A[i].resize( std::min(numNonZero, (uint)_subset.size()) );
+    _B[i].resize( std::min(numNonZero, (uint)_subset.size()) );
+  }
+
+  // generate set from the vector subset
+  set<int> subset_as_set;
+  for (vector<int>::const_iterator i = _subset.begin(); i != _subset.end(); i++ )
+    subset_as_set.insert(*i);
+
+  for (uint dim = 0; dim < this->ui_d; dim++)
+  {
+    double alpha_sum(0.0);
+    double alpha_times_x_sum(0.0);
+
+    uint cntNonzeroFeat(0);
+
+    const multimap< double, SortedVectorSparse<double>::dataelement> & nonzeroElements = this->X_sorted.getFeatureValues(dim).nonzeroElements();
+    // loop through all elements in sorted order
+    for ( SortedVectorSparse<double>::const_elementpointer i = nonzeroElements.begin(); i != nonzeroElements.end(); i++ )
+    {
+      const SortedVectorSparse<double>::dataelement & de = i->second;
+      // index of the feature
+      int index = de.first;
+
+      if (subset_as_set.find(index) != subset_as_set.end()) {
+
+        // transformed element of the feature
+        double elem( de.second );
+
+        alpha_times_x_sum += _alpha[index] * elem;
+        _A[dim][cntNonzeroFeat] = alpha_times_x_sum;
+
+        alpha_sum += _alpha[index];
+        _B[dim][cntNonzeroFeat] = alpha_sum;
+        cntNonzeroFeat++;
+      }
+    }
+  }
+
+}
+
 double *FastMinKernel::hik_prepare_alpha_multiplications_fast(const NICE::VVector & _A,
                                                               const NICE::VVector & _B,
                                                               const Quantization & _q,
@@ -557,6 +610,63 @@ void FastMinKernel::hik_kernel_multiply(const NICE::VVector & _A,
   }
 }
 
+
+void FastMinKernel::hik_kernel_multiply_subset(const NICE::VVector & _A,
+                                        const NICE::VVector & _B,
+                                        const NICE::Vector & _alpha,
+                                        const vector<int> & rowSet,
+                                        const vector<int> & columnSet,
+                                        NICE::Vector & _beta
+                                       ) const
+{
+  _beta.resize( rowSet.size() );
+  _beta.set(0.0);
+
+  map<uint, uint> index_to_subset_index;
+  uint subset_index = 0;
+  for (vector<int>::const_iterator i = rowSet.begin(); i != rowSet.end(); i++, subset_index++)
+    index_to_subset_index.insert( pair<uint, uint>( *i, subset_index ) );
+
+  for (uint dim = 0; dim < this->ui_d; dim++)
+  {
+    const multimap< double, SortedVectorSparse<double>::dataelement> & nonzeroElements = this->X_sorted.getFeatureValues(dim).nonzeroElements();
+
+    uint inversePosition(0);
+    for ( multimap< double, SortedVectorSparse<double>::dataelement>::const_iterator i = nonzeroElements.begin(); i != nonzeroElements.end(); i++)
+    {
+      const SortedVectorSparse<double>::dataelement & de = i->second;
+      uint index = de.first;
+
+      if ( columnSet_as_set.find(index) != columnSet_as_set.end() )
+        inversePosition++;
+
+      map<uint, uint>::const_iterator j = index_to_subset_index.find(index);
+      if (j != index_to_subset_index.end())
+      {
+        uint subset_index = j->second;
+        uint inversePosition = cnt;
+        double fval = de.second;
+
+        double firstPart( _A[dim][inversePosition] );
+
+        // START HERE: tricky thing!
+        double secondPart( _B[dim][this->ui_n-1-nrZeroIndices] - _B[dim][inversePosition]);
+
+        _beta[subset_index] += firstPart + fval * secondPart;
+      }
+    }
+  }
+
+  uint feat = 0;
+  for (map<uint, uint>::const_iterator i = index_to_subset_index.begin(); i != index_to_subset_index.end(); i++)
+  {
+    uint index = i->first;
+    uint subset_index = i->second;
+    _beta[subset_index] += this->d_noise * _alpha[index];
+  }
+}
+
+
 void FastMinKernel::hik_kernel_multiply_fast(const double *_Tlookup,
                                              const Quantization & _q,
                                              const NICE::Vector & _alpha,

FastMinKernel.h

@@ -202,6 +202,17 @@ namespace NICE {
                                              NICE::VVector & _B
                                             ) const;
 
+      /**
+      * @brief  Prepare the efficient HIK-computations part 2 but for only a subset of examples
+      * @author Erik Rodner
+      * @date 15-09-2015 (dd-mm-yyyy)
+      */
+      void hik_prepare_alpha_multiplications_subset(const NICE::Vector & _alpha,
+                                             NICE::VVector & _A,
+                                             NICE::VVector & _B,
+                                             const std::vector<int> & _subset
+                                            ) const;
+
       /**
       * @brief Computing K*alpha with the minimum kernel trick, explicitely exploiting sparsity!!!
       * @author Alexander Freytag

GMHIKernel.cpp

@@ -1,4 +1,4 @@
-/** 
+/**
 * @file GMHIKernel.cpp
 * @brief Fast multiplication with histogram intersection kernel matrices (Implementation)
 * @author Erik Rodner, Alexander Freytag
@@ -39,8 +39,8 @@ void GMHIKernel::multiply (NICE::Vector & y, const NICE::Vector & x) const
     double *T;
     if (useOldPreparation)
     {
-      NICE::VVector A; 
-      NICE::VVector B; 
+      NICE::VVector A;
+      NICE::VVector B;
       // prepare to calculate sum_i x_i K(x,x_i)
       fmk->hik_prepare_alpha_multiplications(x, A, B);
       T = fmk->hik_prepare_alpha_multiplications_fast(A, B, *q, pf);
@@ -49,16 +49,16 @@ void GMHIKernel::multiply (NICE::Vector & y, const NICE::Vector & x) const
     {
       T = fmk->hikPrepareLookupTable(x, *q, pf );
     }
-    fmk->hik_kernel_multiply_fast ( T, *q, x, y ); 
+    fmk->hik_kernel_multiply_fast ( T, *q, x, y );
     delete [] T;
   }
   else //no quantization
   {
-    NICE::VVector A; 
-    NICE::VVector B; 
+    NICE::VVector A;
+    NICE::VVector B;
     // prepare to calculate sum_i x_i K(x,x_i)
     fmk->hik_prepare_alpha_multiplications(x, A, B);
-    
+
     if (verbose)
     {
       int sizeOfDouble (sizeof(double));
@@ -74,7 +74,7 @@ void GMHIKernel::multiply (NICE::Vector & y, const NICE::Vector & x) const
       }
       sizeOfA*=sizeOfDouble;
       sizeOfB*=sizeOfDouble;
-      
+
       std::cerr << "multiplySparse: sizeof(A) + sizeof(B): " << sizeOfA + sizeOfB << std::endl;
     }
     // y = K * x
@@ -84,6 +84,24 @@ void GMHIKernel::multiply (NICE::Vector & y, const NICE::Vector & x) const
   }
 }
 
+/** multiply with a vector: A_subset * x = y */
+void GMHIKernel::multiply (const PartialGenericMatrix::SetType & rowSet, const PartialGenericMatrix::SetType & columnSet, NICE::Vector & y, const NICE::Vector & x) const
+{
+    if (q != NULL) {
+        fthrow(NICE::Exception, "not implemented with quantization")
+    } else {
+      NICE::VVector A;
+      NICE::VVector B;
+      // prepare to calculate sum_i x_i K(x,x_i)
+      fmk->hik_prepare_alpha_multiplications_subset(x, A, B, columnSet);
+
+      // y = K * x
+      //we only need x as input argument to add x*noise to beta
+      //all necessary information for the "real" multiplication is already stored in y
+      fmk->hik_kernel_multiply_subset(A, B, x, rowSet, columnSet, y);
+    }
+}
+
 /** get the number of rows in A */
 uint GMHIKernel::rows () const
 {
@@ -109,7 +127,7 @@ void GMHIKernel::setUseOldPreparation( const bool & _useOldPreparation)
   useOldPreparation = _useOldPreparation;
 }
 
-uint GMHIKernel::getNumParameters() const 
+uint GMHIKernel::getNumParameters() const
 {
   if ( pf == NULL )
     return 0;
@@ -133,7 +151,7 @@ void GMHIKernel::setParameters(const Vector & parameters)
     fthrow(Exception, "Unable to set parameters of a non-parameterized GMHIKernel object");
 
   pf->parameters() = parameters;
-  
+
   fmk->applyFunctionToFeatureMatrix( pf );
 
   // only for debugging with small matrices: fmk->featureMatrix().print();
@@ -146,6 +164,13 @@ void GMHIKernel::getDiagonalElements ( Vector & diagonalElements ) const
   diagonalElements += fmk->getNoise();
 }
 
+double GMHIKernel::getDiagonalElement ( uint i ) const
+{
+    Vector diagonalElements;
+    getDiagonalElements(diagonalElements);
+    return diagonalElements[i];
+}
+
 void GMHIKernel::getFirstDiagonalElement ( double & diagonalElement ) const
 {
   Vector diagonalElements;
@@ -154,12 +179,12 @@ void GMHIKernel::getFirstDiagonalElement ( double & diagonalElement ) const
   // add sigma^2 I
   diagonalElement += fmk->getNoise();
 }
-    
+
 bool GMHIKernel::outOfBounds(const Vector & parameters) const
 {
   if ( pf == NULL && parameters.size() > 0 )
     fthrow(Exception, "Unable to check the bounds of a parameter without any parameterization");
-  
+
   Vector uB = pf->getParameterUpperBounds();
   Vector lB = pf->getParameterLowerBounds();
   if ( uB.size() != parameters.size() || lB.size() != parameters.size() )
@@ -203,8 +228,8 @@ void GMHIKernel::setApproximationScheme(const int & _approxScheme)
 // interface specific methods for incremental extensions
 ///////////////////// INTERFACE ONLINE LEARNABLE /////////////////////
 
-void GMHIKernel::addExample( const NICE::SparseVector * example, 
-			     const double & label, 
+void GMHIKernel::addExample( const NICE::SparseVector * example,
+			     const double & label,
 			     const bool & performOptimizationAfterIncrement
 			   )
 {

GMHIKernel.h

@@ -1,4 +1,4 @@
-/** 
+/**
 * @file GMHIKernel.h
 * @author Erik Rodner, Alexander Freytag
 * @brief Fast multiplication with histogram intersection kernel matrices (Interface)
@@ -11,6 +11,7 @@
 #include <vector>
 
 #include <core/algebra/GenericMatrix.h>
+#include <core/algebra/PartialGenericMatrix.h>
 
 #include "ImplicitKernelMatrix.h"
 #include "FeatureMatrixT.h"
@@ -18,13 +19,13 @@
 
 namespace NICE {
 
- /** 
+ /**
  * @class GMHIKernel
  * @brief Fast multiplication with histogram intersection kernel matrices
  * @author Erik Rodner, Alexander Freytag
  */
 
-class GMHIKernel : public ImplicitKernelMatrix
+class GMHIKernel : public ImplicitKernelMatrix, public PartialGenericMatrix
 {
 
   protected:
@@ -42,10 +43,13 @@ class GMHIKernel : public ImplicitKernelMatrix
 
     /** simple constructor */
     GMHIKernel( FastMinKernel *_fmk, ParameterizedFunction *_pf = NULL, const Quantization *_q = NULL );
-      
+
     /** multiply with a vector: A*x = y */
     virtual void multiply (NICE::Vector & y, const NICE::Vector & x) const;
 
+    /** multiply with a vector: A_subset * x = y */
+    virtual void multiply (const PartialGenericMatrix::SetType & rowSet, const PartialGenericMatrix::SetType & columnSet, NICE::Vector & y, const NICE::Vector & x) const;
+
     /** get the number of rows in A */
     virtual uint rows () const;
 
@@ -54,8 +58,9 @@ class GMHIKernel : public ImplicitKernelMatrix
 
     /** simple destructor */
     virtual ~GMHIKernel();
-   
+
     /** get the diagonal elements of the current matrix */
+    virtual double getDiagonalElement ( uint i ) const;
     virtual void getDiagonalElements ( Vector & diagonalElements ) const;
     virtual void getFirstDiagonalElement ( double & diagonalElement ) const;
 
@@ -69,35 +74,35 @@ class GMHIKernel : public ImplicitKernelMatrix
 
     void setVerbose( const bool & _verbose);
     void setUseOldPreparation( const bool & _useOldPreparation);
-    
+
     virtual double approxFrobNorm() const;
     virtual void setApproximationScheme(const int & _approxScheme);
-    
+
     void setFastMinKernel(NICE::FastMinKernel * _fmk){fmk = _fmk;};
-    
+
     ///////////////////// INTERFACE PERSISTENT /////////////////////
     // interface specific methods for store and restore
     ///////////////////// INTERFACE PERSISTENT /////////////////////
     virtual void restore ( std::istream & is, int format = 0 ) {};//fmk->restore( is, format );};
     virtual void store ( std::ostream & os, int format = 0 ) const {};//fmk->store( os, format );};
     virtual void clear () {};
-    
 
-    
+
+
     ///////////////////// INTERFACE ONLINE LEARNABLE /////////////////////
     // interface specific methods for incremental extensions
-    ///////////////////// INTERFACE ONLINE LEARNABLE /////////////////////    
-    
-    virtual void addExample( const NICE::SparseVector * example, 
-			     const double & label, 
+    ///////////////////// INTERFACE ONLINE LEARNABLE /////////////////////
+
+    virtual void addExample( const NICE::SparseVector * example,
+			     const double & label,
 			     const bool & performOptimizationAfterIncrement = true
 			   );
-			   
+
     virtual void addMultipleExamples( const std::vector< const NICE::SparseVector * > & newExamples,
 				      const NICE::Vector & newLabels,
 				      const bool & performOptimizationAfterIncrement = true
-				    );     
-     
+				    );
+
 };
 
 }