KernelMatrix.h
Go to the documentation of this file.
1 //===========================================================================
2 /*!
3  *
4  *
5  * \brief Kernel Gram matrix
6  *
7  *
8  * \par
9  *
10  *
11  *
12  * \author T. Glasmachers
13  * \date 2007-2012
14  *
15  *
16  * \par Copyright 1995-2015 Shark Development Team
17  *
18  * <BR><HR>
19  * This file is part of Shark.
20  * <http://image.diku.dk/shark/>
21  *
22  * Shark is free software: you can redistribute it and/or modify
23  * it under the terms of the GNU Lesser General Public License as published
24  * by the Free Software Foundation, either version 3 of the License, or
25  * (at your option) any later version.
26  *
27  * Shark is distributed in the hope that it will be useful,
28  * but WITHOUT ANY WARRANTY; without even the implied warranty of
29  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30  * GNU Lesser General Public License for more details.
31  *
32  * You should have received a copy of the GNU Lesser General Public License
33  * along with Shark. If not, see <http://www.gnu.org/licenses/>.
34  *
35  */
36 //===========================================================================
37 
38 
39 #ifndef SHARK_LINALG_KERNELMATRIX_H
40 #define SHARK_LINALG_KERNELMATRIX_H
41 
42 
43 #include <shark/Data/Dataset.h>
44 #include <shark/LinAlg/Base.h>
46 
47 #include <vector>
48 #include <cmath>
49 
50 
51 namespace shark {
52 
53 ///
54 /// \brief Kernel Gram matrix
55 ///
56 /// \par
57 /// The KernelMatrix is the most prominent type of matrix
58 /// for quadratic programming. It provides the Gram matrix
59 /// of a fixed data set with respect to an inner product
60 /// implicitly defined by a kernel function.
61 ///
62 /// \par
63 /// NOTE: The KernelMatrix class stores pointers to the
64 /// data, instead of maintaining a copy of the data. Thus,
65 /// it implicitly assumes that the dataset is not altered
66 /// during the lifetime of the KernelMatrix object. This
67 /// condition is ensured as long as the class is used via
68 /// the various SVM-trainers.
69 ///
70 template <class InputType, class CacheType>
72 {
73 public:
74  typedef CacheType QpFloatType;
75 
76  /// Constructor
77  /// \param kernelfunction kernel function defining the Gram matrix
78  /// \param data data to evaluate the kernel function
80  Data<InputType> const& data)
81  : kernel(kernelfunction)
82  , m_data(data)
83  , m_accessCounter( 0 )
84  {
85  std::size_t elements = m_data.numberOfElements();
86  x.resize(elements);
88  for(std::size_t i = 0; i != elements; ++i,++iter){
89  x[i]=iter.getInnerIterator();
90  }
91  }
92 
93  /// return a single matrix entry
94  QpFloatType operator () (std::size_t i, std::size_t j) const
95  { return entry(i, j); }
96 
97  /// return a single matrix entry
98  QpFloatType entry(std::size_t i, std::size_t j) const
99  {
100  ++m_accessCounter;
101  return (QpFloatType)kernel.eval(*x[i], *x[j]);
102  }
103 
104  /// \brief Computes the i-th row of the kernel matrix.
105  ///
106  ///The entries start,...,end of the i-th row are computed and stored in storage.
107  ///There must be enough room for this operation preallocated.
108  void row(std::size_t i, std::size_t start,std::size_t end, QpFloatType* storage) const{
109  m_accessCounter += end-start;
110 
112  SHARK_PARALLEL_FOR(int j = (int)start; j < (int) end; j++)
113  {
114  storage[j-start] = QpFloatType(kernel.eval(xi, *x[j]));
115  }
116  }
117 
118  /// \brief Computes the kernel-matrix
119  template<class M>
120  void matrix(
122  ) const{
124  }
125 
126  /// swap two variables
127  void flipColumnsAndRows(std::size_t i, std::size_t j){
128  using std::swap;
129  swap(x[i],x[j]);
130  }
131 
132  /// return the size of the quadratic matrix
133  std::size_t size() const
134  { return x.size(); }
135 
136  /// query the kernel access counter
137  unsigned long long getAccessCount() const
138  { return m_accessCounter; }
139 
140  /// reset the kernel access counter
142  { m_accessCounter = 0; }
143 
144 protected:
145  /// Kernel function defining the kernel Gram matrix
147 
149 
151  /// Array of data pointers for kernel evaluations
152  std::vector<PointerType> x;
153 
154  /// counter for the kernel accesses
155  mutable unsigned long long m_accessCounter;
156 };
157 
158 //~ ///\brief Specialization for dense vectors which often can be computed much faster
159 //~ template <class T, class CacheType>
160 //~ class KernelMatrix<blas::vector<T>, CacheType>
161 //~ {
162 //~ public:
163 
164  //~ //////////////////////////////////////////////////////////////////
165  //~ // The types below define the type used for caching kernel values. The default is float,
166  //~ // since this type offers sufficient accuracy in the vast majority of cases, at a memory
167  //~ // cost of only four bytes. However, the type definition makes it easy to use double instead
168  //~ // (e.g., in case high accuracy training is needed).
169  //~ typedef CacheType QpFloatType;
170  //~ typedef blas::vector<T> InputType;
171 
172  //~ /// Constructor
173  //~ /// \param kernelfunction kernel function defining the Gram matrix
174  //~ /// \param data data to evaluate the kernel function
175  //~ KernelMatrix(
176  //~ AbstractKernelFunction<InputType> const& kernelfunction,
177  //~ Data<InputType> const& data)
178  //~ : kernel(kernelfunction)
179  //~ , m_data(data)
180  //~ , m_batchStart(data.numberOfBatches())
181  //~ , m_accessCounter( 0 )
182  //~ {
183  //~ m_data.makeIndependent();
184  //~ std::size_t elements = m_data.numberOfElements();
185  //~ x.resize(elements);
186  //~ typename Data<InputType>::element_range::iterator iter=m_data.elements().begin();
187  //~ for(std::size_t i = 0; i != elements; ++i,++iter){
188  //~ x[i]=iter.getInnerIterator();
189  //~ }
190 
191  //~ for(std::size_t i = 0,start = 0; i != m_data.numberOfBatches(); ++i){
192  //~ m_batchStart[i] = start;
193  //~ start+= m_data.batch(i).size1();
194  //~ }
195  //~ }
196 
197  //~ /// return a single matrix entry
198  //~ QpFloatType operator () (std::size_t i, std::size_t j) const
199  //~ { return entry(i, j); }
200 
201  //~ /// return a single matrix entry
202  //~ QpFloatType entry(std::size_t i, std::size_t j) const
203  //~ {
204  //~ ++m_accessCounter;
205  //~ return (QpFloatType)kernel.eval(*x[i], *x[j]);
206  //~ }
207 
208  //~ /// \brief Computes the i-th row of the kernel matrix.
209  //~ ///
210  //~ ///The entries start,...,end of the i-th row are computed and stored in storage.
211  //~ ///There must be enough room for this operation preallocated.
212  //~ void row(std::size_t k, std::size_t start,std::size_t end, QpFloatType* storage) const
213  //~ {
214  //~ m_accessCounter +=end-start;
215 
216  //~ typename AbstractKernelFunction<InputType>::ConstInputReference xi = *x[k];
217  //~ typename blas::matrix<T> mx(1,xi.size());
218  //~ noalias(blas::row(mx,0))=xi;
219 
220  //~ int numBatches = (int)m_data.numberOfBatches();
221  //~ SHARK_PARALLEL_FOR(int i = 0; i < numBatches; i++)
222  //~ {
223  //~ std::size_t pos = m_batchStart[i];
224  //~ std::size_t batchSize = m_data.batch(i).size1();
225  //~ if(!(pos+batchSize < start || pos > end)){
226  //~ RealMatrix rowpart(1,batchSize);
227  //~ kernel.eval(mx,m_data.batch(i),rowpart);
228  //~ std::size_t batchStart = (start <=pos) ? 0: start-pos;
229  //~ std::size_t batchEnd = (pos+batchSize > end) ? end-pos: batchSize;
230  //~ for(std::size_t j = batchStart; j != batchEnd;++j){
231  //~ storage[pos+j-start] = static_cast<QpFloatType>(rowpart(0,j));
232  //~ }
233  //~ }
234  //~ }
235  //~ }
236 
237  //~ /// \brief Computes the kernel-matrix
238  //~ template<class M>
239  //~ void matrix(
240  //~ blas::matrix_expression<M> & storage
241  //~ ) const{
242  //~ calculateRegularizedKernelMatrix(kernel,m_data,storage);
243  //~ }
244 
245  //~ /// swap two variables
246  //~ void flipColumnsAndRows(std::size_t i, std::size_t j){
247  //~ if( i == j ) return;
248  //~ swap(*x[i],*x[j]);
249  //~ }
250 
251  //~ /// return the size of the quadratic matrix
252  //~ std::size_t size() const
253  //~ { return x.size(); }
254 
255  //~ /// query the kernel access counter
256  //~ unsigned long long getAccessCount() const
257  //~ { return m_accessCounter; }
258 
259  //~ /// reset the kernel access counter
260  //~ void resetAccessCount()
261  //~ { m_accessCounter = 0; }
262 
263 //~ protected:
264  //~ /// Kernel function defining the kernel Gram matrix
265  //~ const AbstractKernelFunction<InputType>& kernel;
266 
267  //~ Data<InputType> m_data;
268 
269  //~ typedef typename Batch<InputType>::iterator PointerType;
270  //~ /// Array of data pointers for kernel evaluations
271  //~ std::vector<PointerType> x;
272 
273  //~ std::vector<std::size_t> m_batchStart;
274 
275  //~ mutable unsigned long long m_accessCounter;
276 //~ };
277 
278 }
279 #endif