Tapkee
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
methods.hpp
Go to the documentation of this file.
1 /* This software is distributed under BSD 3-clause license (see LICENSE file).
2  *
3  * Copyright (c) 2012-2013 Sergey Lisitsyn, Fernando Iglesias
4  */
5 
6 #ifndef TAPKEE_METHODS_H_
7 #define TAPKEE_METHODS_H_
8 
9 /* Tapkee includes */
10 #include <tapkee/defines.hpp>
11 #include <tapkee/utils/naming.hpp>
12 #include <tapkee/utils/time.hpp>
13 #include <tapkee/utils/logging.hpp>
25 #include <tapkee/routines/pca.hpp>
27 #include <tapkee/routines/spe.hpp>
28 #include <tapkee/routines/fa.hpp>
32 /* End of Tapkee includes */
33 
34 namespace tapkee
35 {
37 namespace tapkee_internal
38 {
39 
40 template <class RandomAccessIterator, class KernelCallback,
41  class DistanceCallback, class FeaturesCallback>
43 {
44 public:
45 
46  ImplementationBase(RandomAccessIterator b, RandomAccessIterator e,
47  KernelCallback k, DistanceCallback d, FeaturesCallback f,
48  ParametersSet& pmap, const Context& ctx) :
49  parameters(pmap), context(ctx), kernel(k), distance(d), features(f),
50  plain_distance(PlainDistance<RandomAccessIterator,DistanceCallback>(distance)),
51  kernel_distance(KernelDistance<RandomAccessIterator,KernelCallback>(kernel)),
52  begin(b), end(e),
58  {
59  n_vectors = (end-begin);
60 
61  target_dimension = parameters(keywords::target_dimension);
62  n_neighbors = parameters(keywords::num_neighbors).checked().positive();
63 
64  if (n_vectors > 0)
65  {
67  .inRange(static_cast<IndexType>(1),static_cast<IndexType>(n_vectors));
69  .inRange(static_cast<IndexType>(3),static_cast<IndexType>(n_vectors));
70  }
71 
72  eigen_method = parameters(keywords::eigen_method);
73  neighbors_method = parameters(keywords::neighbors_method);
74  check_connectivity = parameters(keywords::check_connectivity);
75  width = parameters(keywords::gaussian_kernel_width).checked().positive();
76  timesteps = parameters(keywords::diffusion_map_timesteps).checked().positive();
77  eigenshift = parameters(keywords::nullspace_shift);
78  traceshift = parameters(keywords::klle_shift);
79  max_iteration = parameters(keywords::max_iteration);
80  tolerance = parameters(keywords::spe_tolerance).checked().positive();
81  n_updates = parameters(keywords::spe_num_updates).checked().positive();
82  theta = parameters(keywords::sne_theta).checked().nonNegative();
83  squishing_rate = parameters(keywords::squishing_rate);
84  global_strategy = parameters(keywords::spe_global_strategy);
85  epsilon = parameters(keywords::fa_epsilon).checked().nonNegative();
86  perplexity = parameters(keywords::sne_perplexity).checked().nonNegative();
87  ratio = parameters(keywords::landmark_ratio);
88 
90  {
91  current_dimension = features.dimension();
92  }
93  else
94  {
96  }
97  }
98 
100  {
101  if (context.is_cancelled())
102  throw cancelled_exception();
103 
104  using std::mem_fun_ref_t;
105  using std::mem_fun_ref;
106  typedef std::mem_fun_ref_t<TapkeeOutput,ImplementationBase> ImplRef;
107 
108 #define tapkee_method_handle(X) \
109  case X: \
110  { \
111  timed_context tctx__("[+] embedding with " # X); \
112  ImplRef ref = conditional_select< \
113  ((!MethodTraits<X>::needs_kernel) || (!is_dummy<KernelCallback>::value)) && \
114  ((!MethodTraits<X>::needs_distance) || (!is_dummy<DistanceCallback>::value)) && \
115  ((!MethodTraits<X>::needs_features) || (!is_dummy<FeaturesCallback>::value)), \
116  ImplRef>()(mem_fun_ref(&ImplementationBase::embed##X), \
117  mem_fun_ref(&ImplementationBase::embedEmpty)); \
118  return ref(*this); \
119  } \
120  break \
121 
122  switch (method)
123  {
144  }
145 #undef tapkee_method_handle
146  return TapkeeOutput();
147  }
148 
149 private:
150 
151  static const IndexType SkipOneEigenvalue = 1;
152  static const IndexType SkipNoEigenvalues = 0;
153 
156  KernelCallback kernel;
157  DistanceCallback distance;
158  FeaturesCallback features;
161 
162  RandomAccessIterator begin;
163  RandomAccessIterator end;
164 
183 
186 
187  template<class Distance>
189  {
191  }
192 
194  {
196  }
197 
199  {
200  throw unsupported_method_error("Some callback is missed");
201  return TapkeeOutput();
202  }
203 
205  {
207  SparseWeightMatrix weight_matrix =
209  DenseMatrix embedding =
210  eigendecomposition<SparseWeightMatrix,SparseInverseMatrixOperation>(eigen_method,
211  weight_matrix,target_dimension,SkipOneEigenvalue).first;
212 
213  return TapkeeOutput(embedding, unimplementedProjectingFunction());
214  }
215 
217  {
219  SparseWeightMatrix weight_matrix =
221  DenseMatrix embedding =
222  eigendecomposition<SparseWeightMatrix,SparseInverseMatrixOperation>(eigen_method,
223  weight_matrix,target_dimension,SkipOneEigenvalue).first;
224 
225  return TapkeeOutput(embedding, unimplementedProjectingFunction());
226  }
227 
229  {
230  #ifdef TAPKEE_GPU
231  #define DM_MATRIX_OP GPUDenseImplicitSquareMatrixOperation
232  #else
233  #define DM_MATRIX_OP DenseImplicitSquareSymmetricMatrixOperation
234  #endif
235 
236  DenseSymmetricMatrix diffusion_matrix =
238  DenseMatrix embedding =
239  eigendecomposition<DenseSymmetricMatrix,DM_MATRIX_OP>(eigen_method,diffusion_matrix,
241 
242  return TapkeeOutput(embedding, unimplementedProjectingFunction());
243 
244  #undef DM_MATRIX_OP
245  }
246 
248  {
249  #ifdef TAPKEE_GPU
250  #define MDS_MATRIX_OP GPUDenseImplicitSquareMatrixOperation
251  #else
252  #define MDS_MATRIX_OP DenseMatrixOperation
253  #endif
254 
256  centerMatrix(distance_matrix);
257  distance_matrix.array() *= -0.5;
258  EigendecompositionResult embedding =
259  eigendecomposition<DenseSymmetricMatrix,MDS_MATRIX_OP>(eigen_method,
260  distance_matrix,target_dimension,SkipNoEigenvalues);
261 
262  for (IndexType i=0; i<static_cast<IndexType>(target_dimension); i++)
263  embedding.first.col(i).array() *= sqrt(embedding.second(i));
264  return TapkeeOutput(embedding.first, unimplementedProjectingFunction());
265  #undef MDS_MATRIX_OP
266  }
267 
269  {
270  ratio.checked()
271  .inClosedRange(static_cast<ScalarType>(3.0/n_vectors),
272  static_cast<ScalarType>(1.0));
273 
274  Landmarks landmarks =
276  DenseSymmetricMatrix distance_matrix =
278  DenseVector landmark_distances_squared = distance_matrix.colwise().mean();
279  centerMatrix(distance_matrix);
280  distance_matrix.array() *= -0.5;
281  EigendecompositionResult landmarks_embedding =
282  eigendecomposition<DenseSymmetricMatrix,DenseMatrixOperation>(eigen_method,
283  distance_matrix,target_dimension,SkipNoEigenvalues);
284  for (IndexType i=0; i<static_cast<IndexType>(target_dimension); i++)
285  landmarks_embedding.first.col(i).array() *= sqrt(landmarks_embedding.second(i));
286  return TapkeeOutput(triangulate(begin,end,distance,landmarks,
287  landmark_distances_squared,landmarks_embedding,target_dimension), unimplementedProjectingFunction());
288  }
289 
291  {
293  DenseSymmetricMatrix shortest_distances_matrix =
295  shortest_distances_matrix = shortest_distances_matrix.array().square();
296  centerMatrix(shortest_distances_matrix);
297  shortest_distances_matrix.array() *= -0.5;
298 
299  EigendecompositionResult embedding =
300  eigendecomposition<DenseSymmetricMatrix,DenseMatrixOperation>(eigen_method,
301  shortest_distances_matrix,target_dimension,SkipNoEigenvalues);
302 
303  for (IndexType i=0; i<static_cast<IndexType>(target_dimension); i++)
304  embedding.first.col(i).array() *= sqrt(embedding.second(i));
305 
306  return TapkeeOutput(embedding.first, unimplementedProjectingFunction());
307  }
308 
310  {
311  ratio.checked()
312  .inClosedRange(static_cast<ScalarType>(3.0/n_vectors),
313  static_cast<ScalarType>(1.0));
314 
316  Landmarks landmarks =
318  DenseMatrix distance_matrix =
320  distance_matrix = distance_matrix.array().square();
321 
322  DenseVector col_means = distance_matrix.colwise().mean();
323  DenseVector row_means = distance_matrix.rowwise().mean();
324  ScalarType grand_mean = distance_matrix.mean();
325  distance_matrix.array() += grand_mean;
326  distance_matrix.colwise() -= row_means;
327  distance_matrix.rowwise() -= col_means.transpose();
328  distance_matrix.array() *= -0.5;
329 
330  EigendecompositionResult landmarks_embedding;
331 
332  if (eigen_method.is(Dense))
333  {
334  DenseMatrix distance_matrix_sym = distance_matrix*distance_matrix.transpose();
335  landmarks_embedding = eigendecomposition<DenseSymmetricMatrix,DenseImplicitSquareMatrixOperation>
336  (eigen_method,distance_matrix_sym,target_dimension,SkipNoEigenvalues);
337  }
338  else
339  {
340  landmarks_embedding = eigendecomposition<DenseSymmetricMatrix,DenseImplicitSquareMatrixOperation>
342  }
343 
344  DenseMatrix embedding = distance_matrix.transpose()*landmarks_embedding.first;
345 
346  for (IndexType i=0; i<static_cast<IndexType>(target_dimension); i++)
347  embedding.col(i).array() /= sqrt(sqrt(landmarks_embedding.second(i)));
348  return TapkeeOutput(embedding,unimplementedProjectingFunction());
349  }
350 
352  {
354  SparseWeightMatrix weight_matrix =
356  DenseSymmetricMatrixPair eig_matrices =
359  EigendecompositionResult projection_result =
360  generalized_eigendecomposition<DenseSymmetricMatrix,DenseSymmetricMatrix,DenseInverseMatrixOperation>(
361  eigen_method,eig_matrices.first,eig_matrices.second,target_dimension,SkipNoEigenvalues);
362  DenseVector mean_vector =
364  tapkee::ProjectingFunction projecting_function(new tapkee::MatrixProjectionImplementation(projection_result.first,mean_vector));
365  return TapkeeOutput(project(projection_result.first,mean_vector,begin,end,features,current_dimension),projecting_function);
366  }
367 
369  {
371  SparseWeightMatrix weight_matrix =
373  return TapkeeOutput(eigendecomposition<SparseWeightMatrix,SparseInverseMatrixOperation>(eigen_method,
375  }
376 
378  {
380  Laplacian laplacian =
382  return TapkeeOutput(generalized_eigendecomposition<SparseWeightMatrix,DenseDiagonalMatrix,SparseInverseMatrixOperation>(
383  eigen_method,laplacian.first,laplacian.second,target_dimension,SkipOneEigenvalue).first, unimplementedProjectingFunction());
384  }
385 
387  {
389  Laplacian laplacian =
391  DenseSymmetricMatrixPair eigenproblem_matrices =
392  construct_locality_preserving_eigenproblem(laplacian.first,laplacian.second,begin,end,
394  EigendecompositionResult projection_result =
395  generalized_eigendecomposition<DenseSymmetricMatrix,DenseSymmetricMatrix,DenseInverseMatrixOperation>(
396  eigen_method,eigenproblem_matrices.first,eigenproblem_matrices.second,target_dimension,SkipNoEigenvalues);
397  DenseVector mean_vector =
399  tapkee::ProjectingFunction projecting_function(new tapkee::MatrixProjectionImplementation(projection_result.first,mean_vector));
400  return TapkeeOutput(project(projection_result.first,mean_vector,begin,end,features,current_dimension), projecting_function);
401  }
402 
404  {
405  DenseVector mean_vector =
407  DenseSymmetricMatrix centered_covariance_matrix =
409  EigendecompositionResult projection_result =
410  eigendecomposition<DenseSymmetricMatrix,DenseMatrixOperation>(eigen_method,centered_covariance_matrix,target_dimension,SkipNoEigenvalues);
411  tapkee::ProjectingFunction projecting_function(new tapkee::MatrixProjectionImplementation(projection_result.first,mean_vector));
412  return TapkeeOutput(project(projection_result.first,mean_vector,begin,end,features,current_dimension), projecting_function);
413  }
414 
416  {
417  DenseMatrix projection_matrix =
419  DenseVector mean_vector =
421 
422  tapkee::ProjectingFunction projecting_function(new tapkee::MatrixProjectionImplementation(projection_matrix,mean_vector));
423  return TapkeeOutput(project(projection_matrix,mean_vector,begin,end,features,current_dimension), projecting_function);
424  }
425 
427  {
428  DenseSymmetricMatrix centered_kernel_matrix =
430  EigendecompositionResult embedding = eigendecomposition<DenseSymmetricMatrix,DenseMatrixOperation>(eigen_method,
431  centered_kernel_matrix,target_dimension,SkipNoEigenvalues);
432  for (IndexType i=0; i<static_cast<IndexType>(target_dimension); i++)
433  embedding.first.col(i).array() *= sqrt(embedding.second(i));
434  return TapkeeOutput(embedding.first, unimplementedProjectingFunction());
435  }
436 
438  {
440  SparseWeightMatrix weight_matrix =
442  DenseSymmetricMatrixPair eig_matrices =
443  construct_lltsa_eigenproblem(weight_matrix,begin,end,
445  EigendecompositionResult projection_result =
446  generalized_eigendecomposition<DenseSymmetricMatrix,DenseSymmetricMatrix,DenseInverseMatrixOperation>(
447  eigen_method,eig_matrices.first,eig_matrices.second,target_dimension,SkipNoEigenvalues);
448  DenseVector mean_vector =
450  tapkee::ProjectingFunction projecting_function(new tapkee::MatrixProjectionImplementation(projection_result.first,mean_vector));
451  return TapkeeOutput(project(projection_result.first,mean_vector,begin,end,features,current_dimension),
452  projecting_function);
453  }
454 
456  {
457  Neighbors neighbors;
458  if (global_strategy.is(false))
459  {
460  neighbors = findNeighborsWith(plain_distance);
461  }
462 
463  return TapkeeOutput(spe_embedding(begin,end,distance,neighbors,
465  }
466 
468  {
469  DenseMatrix feature_matrix =
471  return TapkeeOutput(feature_matrix.transpose(),tapkee::ProjectingFunction());
472  }
473 
475  {
479  }
480 
482  {
484  .inClosedRange(static_cast<ScalarType>(0.0),
485  static_cast<ScalarType>((n_vectors-1)/3.0));
486 
487  DenseMatrix data =
489 
490  DenseMatrix embedding(static_cast<IndexType>(target_dimension),n_vectors);
491  tsne::TSNE tsne;
492  tsne.run(data.data(),n_vectors,current_dimension,embedding.data(),target_dimension,perplexity,theta);
493 
494  return TapkeeOutput(embedding.transpose(), unimplementedProjectingFunction());
495  }
496 
498  {
500  .inRange(static_cast<ScalarType>(0.0),
501  static_cast<ScalarType>(1.0));
502 
503  DenseMatrix embedding =
505 
507 
509 
510  return TapkeeOutput(embedding, tapkee::ProjectingFunction());
511  }
512 
513 };
514 
515 template <class RandomAccessIterator, class KernelCallback,
516  class DistanceCallback, class FeaturesCallback>
517 ImplementationBase<RandomAccessIterator,KernelCallback,DistanceCallback,FeaturesCallback>
518  initialize(RandomAccessIterator begin, RandomAccessIterator end,
519  KernelCallback kernel, DistanceCallback distance, FeaturesCallback features,
520  ParametersSet& pmap, const Context& ctx)
521 {
523  begin,end,kernel,distance,features,pmap,ctx);
524 }
525 
526 } // End of namespace tapkee_internal
527 } // End of namespace tapkee
528 
529 #endif
ImplementationBase(RandomAccessIterator b, RandomAccessIterator e, KernelCallback k, DistanceCallback d, FeaturesCallback f, ParametersSet &pmap, const Context &ctx)
Definition: methods.hpp:46
void run(double *X, int N, int D, double *Y, int no_dims, double perplexity, double theta)
Definition: tsne.hpp:58
ScalarType distance(Callback &cb, const CoverTreePoint< RandomAccessIterator > &l, const CoverTreePoint< RandomAccessIterator > &r, ScalarType upper_bound)
void centerMatrix(DenseMatrix &matrix)
Definition: matrix.hpp:9
Neighbors find_neighbors(NeighborsMethod method, const RandomAccessIterator &begin, const RandomAccessIterator &end, const Callback &callback, IndexType k, bool check_connectivity)
Definition: neighbors.hpp:173
Eigen::Matrix< tapkee::ScalarType, Eigen::Dynamic, Eigen::Dynamic > DenseMatrix
dense matrix type (non-overridable)
Definition: types.hpp:23
DenseSymmetricMatrix compute_centered_kernel_matrix(RandomAccessIterator begin, RandomAccessIterator end, KernelCallback callback)
Definition: pca.hpp:76
void manifold_sculpting_embed(RandomAccessIterator begin, RandomAccessIterator end, DenseMatrix &data, IndexType target_dimension, const Neighbors &neighbors, DistanceCallback callback, IndexType max_iteration, ScalarType squishing_rate)
DenseSymmetricMatrix compute_covariance_matrix(RandomAccessIterator begin, RandomAccessIterator end, const DenseVector &mean, FeatureVectorCallback callback, IndexType dimension)
Definition: pca.hpp:57
DenseSymmetricMatrixPair construct_lltsa_eigenproblem(SparseWeightMatrix W, RandomAccessIterator begin, RandomAccessIterator end, FeatureVectorCallback feature_vector_callback, IndexType dimension)
static const IndexType SkipOneEigenvalue
Definition: methods.hpp:151
SparseWeightMatrix hessian_weight_matrix(RandomAccessIterator begin, RandomAccessIterator end, const Neighbors &neighbors, PairwiseCallback callback, const IndexType target_dimension)
TAPKEE_INTERNAL_PAIR< tapkee::SparseWeightMatrix, tapkee::DenseDiagonalMatrix > Laplacian
Definition: synonyms.hpp:43
Eigen library dense method (could be useful for debugging). Computes all eigenvectors thus can be ver...
static const IndexType SkipNoEigenvalues
Definition: methods.hpp:152
DenseSymmetricMatrix compute_distance_matrix(RandomAccessIterator begin, RandomAccessIterator, const Landmarks &landmarks, PairwiseCallback callback)
CheckedParameter checked()
Definition: parameters.hpp:300
double ScalarType
default scalar value (can be overrided with TAPKEE_CUSTOM_INTERNAL_NUMTYPE define) ...
Definition: types.hpp:15
TAPKEE_INTERNAL_PAIR< tapkee::DenseSymmetricMatrix, tapkee::DenseSymmetricMatrix > DenseSymmetricMatrixPair
Definition: synonyms.hpp:44
DenseSymmetricMatrixPair construct_neighborhood_preserving_eigenproblem(SparseWeightMatrix W, RandomAccessIterator begin, RandomAccessIterator end, FeatureVectorCallback feature_vector_callback, IndexType dimension)
An exception type that is thrown when computations were cancelled.
Definition: exceptions.hpp:76
Eigen::SparseMatrix< tapkee::ScalarType > SparseWeightMatrix
sparse weight matrix type (non-overridable)
Definition: types.hpp:29
DenseSymmetricMatrix compute_diffusion_matrix(RandomAccessIterator begin, RandomAccessIterator end, DistanceCallback callback, const IndexType timesteps, const ScalarType width)
Computes diffusion process matrix. Uses the following algorithm:
TapkeeOutput embedtDistributedStochasticNeighborEmbedding()
Definition: methods.hpp:481
DenseSymmetricMatrix compute_shortest_distances_matrix(const RandomAccessIterator &begin, const RandomAccessIterator &end, const Neighbors &neighbors, DistanceCallback callback)
Computes shortest distances (so-called geodesic distances) using Dijkstra algorithm.
Definition: isomap.hpp:44
Eigen::Matrix< tapkee::ScalarType, Eigen::Dynamic, 1 > DenseVector
dense vector type (non-overridable)
Definition: types.hpp:21
DenseMatrix triangulate(RandomAccessIterator begin, RandomAccessIterator end, PairwiseCallback distance_callback, const Landmarks &landmarks, const DenseVector &landmark_distances_squared, EigendecompositionResult &landmarks_embedding, IndexType target_dimension)
DenseMatrix dense_matrix_from_features(FeaturesCallback features, IndexType dimension, RandomAccessIterator begin, RandomAccessIterator end)
Definition: features.hpp:19
Laplacian compute_laplacian(RandomAccessIterator begin, RandomAccessIterator end, const Neighbors &neighbors, DistanceCallback callback, ScalarType width)
Computes laplacian of neighborhood graph.
int IndexType
indexing type (non-overridable) set to int for compatibility with OpenMP 2.0
Definition: types.hpp:19
TapkeeOutput embedUsing(DimensionReductionMethod method)
Definition: methods.hpp:99
PlainDistance< RandomAccessIterator, DistanceCallback > plain_distance
Definition: methods.hpp:159
DimensionReductionMethod
Dimension reduction methods.
DenseVector compute_mean(RandomAccessIterator begin, RandomAccessIterator end, FeatureVectorCallback callback, IndexType dimension)
Definition: pca.hpp:42
DenseMatrix project(RandomAccessIterator begin, RandomAccessIterator end, FeatureVectorCallback callback, IndexType dimension, const IndexType max_iter, const ScalarType epsilon, const IndexType target_dimension, const DenseVector &mean_vector)
Definition: fa.hpp:15
A pimpl wrapper for projecting function.
Definition: projection.hpp:25
static tapkee::ProjectingFunction unimplementedProjectingFunction()
Definition: methods.hpp:193
Return result of the library - a pair of DenseMatrix (embedding) and ProjectingFunction.
Definition: defines.hpp:43
CheckedParameter & inRange(T lower, T upper)
Definition: parameters.hpp:244
TAPKEE_INTERNAL_VECTOR< tapkee::tapkee_internal::LocalNeighbors > Neighbors
Definition: synonyms.hpp:40
SparseWeightMatrix linear_weight_matrix(const RandomAccessIterator &begin, const RandomAccessIterator &end, const Neighbors &neighbors, PairwiseCallback callback, const ScalarType shift, const ScalarType trace_shift)
#define tapkee_method_handle(X)
ImplementationBase< RandomAccessIterator, KernelCallback, DistanceCallback, FeaturesCallback > initialize(RandomAccessIterator begin, RandomAccessIterator end, KernelCallback kernel, DistanceCallback distance, FeaturesCallback features, ParametersSet &pmap, const Context &ctx)
Definition: methods.hpp:518
An exception type that is thrown when unsupported method is called.
Definition: exceptions.hpp:47
TAPKEE_INTERNAL_PAIR< tapkee::DenseMatrix, tapkee::DenseVector > EigendecompositionResult
Definition: synonyms.hpp:41
Landmarks select_landmarks_random(RandomAccessIterator begin, RandomAccessIterator end, ScalarType ratio)
TAPKEE_INTERNAL_VECTOR< tapkee::IndexType > Landmarks
Definition: synonyms.hpp:42
KernelDistance< RandomAccessIterator, KernelCallback > kernel_distance
Definition: methods.hpp:160
tapkee::DenseMatrix DenseSymmetricMatrix
dense symmetric matrix (non-overridable, currently just dense matrix, can be improved later) ...
Definition: types.hpp:25
DenseSymmetricMatrixPair construct_locality_preserving_eigenproblem(SparseWeightMatrix &L, DenseDiagonalMatrix &D, RandomAccessIterator begin, RandomAccessIterator end, FeatureVectorCallback feature_vector_callback, IndexType dimension)
CheckedParameter & inClosedRange(T lower, T upper)
Definition: parameters.hpp:258
DenseMatrix gaussian_projection_matrix(IndexType target_dimension, IndexType current_dimension)
DenseMatrix spe_embedding(RandomAccessIterator begin, RandomAccessIterator end, PairwiseCallback callback, const Neighbors &neighbors, IndexType target_dimension, bool global_strategy, ScalarType tolerance, int nupdates, IndexType max_iter)
Definition: spe.hpp:20
Basic ProjectionImplementation that subtracts mean from the vector and multiplies projecting matrix w...
Definition: projection.hpp:46
SparseWeightMatrix tangent_weight_matrix(RandomAccessIterator begin, RandomAccessIterator end, const Neighbors &neighbors, PairwiseCallback callback, const IndexType target_dimension, const ScalarType shift, const bool partial_eigendecomposer=false)