6 #ifndef TAPKEE_METHODS_H_ 7 #define TAPKEE_METHODS_H_ 37 namespace tapkee_internal
69 InRange(T l, T u) : lower(l), upper(u) { }
72 return (v>=lower) && (v<upper);
88 return (v>=lower) && (v<=upper);
98 template <
class RandomAccessIterator,
class KernelCallback,
99 class DistanceCallback,
class FeaturesCallback>
105 KernelCallback k, DistanceCallback d, FeaturesCallback f,
107 parameters(pmap), context(ctx), kernel(k),
distance(d), features(f),
109 kernel_distance(
KernelDistance<RandomAccessIterator,KernelCallback>(kernel)),
110 begin(b), end(e), p_computation_strategy(),
111 p_eigen_method(), p_neighbors_method(), p_eigenshift(), p_traceshift(),
112 p_check_connectivity(), p_n_neighbors(), p_width(), p_timesteps(),
113 p_ratio(), p_max_iteration(), p_tolerance(), p_n_updates(), p_perplexity(),
114 p_theta(), p_squishing_rate(), p_global_strategy(), p_epsilon(), p_target_dimension(),
115 n_vectors(0), current_dimension(0)
117 n_vectors = (end-begin);
119 p_target_dimension = parameters[target_dimension];
131 p_computation_strategy = parameters[computation_strategy];
132 p_eigen_method = parameters[eigen_method];
133 p_neighbors_method = parameters[neighbors_method];
134 p_check_connectivity = parameters[check_connectivity];
137 p_eigenshift = parameters[nullspace_shift];
138 p_traceshift = parameters[klle_shift];
139 p_max_iteration = parameters[max_iteration];
143 p_squishing_rate = parameters[squishing_rate];
144 p_global_strategy = parameters[spe_global_strategy];
147 p_ratio = parameters[landmark_ratio];
150 current_dimension = features.dimension();
152 current_dimension = 0;
157 if (context.is_cancelled())
160 using std::mem_fun_ref_t;
161 using std::mem_fun_ref;
162 typedef std::mem_fun_ref_t<TapkeeOutput,ImplementationBase> ImplRef;
164 #define tapkee_method_handle(X) \ 167 timed_context tctx__("[+] embedding with " # X); \ 168 ImplRef ref = conditional_select< \ 169 ((!MethodTraits<X>::needs_kernel) || (!is_dummy<KernelCallback>::value)) && \ 170 ((!MethodTraits<X>::needs_distance) || (!is_dummy<DistanceCallback>::value)) && \ 171 ((!MethodTraits<X>::needs_features) || (!is_dummy<FeaturesCallback>::value)), \ 172 ImplRef>()(mem_fun_ref(&ImplementationBase::embed##X), \ 173 mem_fun_ref(&ImplementationBase::embedEmpty)); \ 201 #undef tapkee_method_handle 241 template<
class Distance>
244 return find_neighbors(p_neighbors_method,begin,end,d,p_n_neighbors,p_check_connectivity);
260 Neighbors neighbors = findNeighborsWith(kernel_distance);
265 weight_matrix,p_target_dimension).first;
267 return TapkeeOutput(embedding, unimplementedProjectingFunction());
272 Neighbors neighbors = findNeighborsWith(kernel_distance);
277 weight_matrix,p_target_dimension).first;
279 return TapkeeOutput(embedding, unimplementedProjectingFunction());
290 DenseMatrix embedding = (decomposition_result.first).leftCols(target_dimension);
292 for (
IndexType i=0; i<target_dimension; i++)
293 embedding.col(i).array() *= pow(decomposition_result.second(i),
static_cast<IndexType>(p_timesteps));
295 for (
IndexType i=0; i<target_dimension; i++)
296 embedding.col(i).array() /= decomposition_result.first.col(target_dimension).array();
297 return TapkeeOutput(embedding, unimplementedProjectingFunction());
304 distance_matrix.array() *= -0.5;
307 distance_matrix,p_target_dimension);
309 for (
IndexType i=0; i<static_cast<IndexType>(p_target_dimension); i++)
310 embedding.first.col(i).array() *= sqrt(embedding.second(i));
311 return TapkeeOutput(embedding.first, unimplementedProjectingFunction());
323 DenseVector landmark_distances_squared = distance_matrix.colwise().mean();
325 distance_matrix.array() *= -0.5;
328 distance_matrix,p_target_dimension);
329 for (
IndexType i=0; i<static_cast<IndexType>(p_target_dimension); i++)
330 landmarks_embedding.first.col(i).array() *= sqrt(landmarks_embedding.second(i));
332 landmark_distances_squared,landmarks_embedding,p_target_dimension), unimplementedProjectingFunction());
337 Neighbors neighbors = findNeighborsWith(plain_distance);
340 shortest_distances_matrix = shortest_distances_matrix.array().square();
342 shortest_distances_matrix.array() *= -0.5;
346 shortest_distances_matrix,p_target_dimension);
348 for (
IndexType i=0; i<static_cast<IndexType>(p_target_dimension); i++)
349 embedding.first.col(i).array() *= sqrt(embedding.second(i));
351 return TapkeeOutput(embedding.first, unimplementedProjectingFunction());
358 Neighbors neighbors = findNeighborsWith(plain_distance);
363 distance_matrix = distance_matrix.array().square();
365 DenseVector col_means = distance_matrix.colwise().mean();
366 DenseVector row_means = distance_matrix.rowwise().mean();
367 ScalarType grand_mean = distance_matrix.mean();
368 distance_matrix.array() += grand_mean;
369 distance_matrix.colwise() -= row_means;
370 distance_matrix.rowwise() -= col_means.transpose();
371 distance_matrix.array() *= -0.5;
377 DenseMatrix distance_matrix_sym = distance_matrix*distance_matrix.transpose();
387 DenseMatrix embedding = distance_matrix.transpose()*landmarks_embedding.first;
389 for (
IndexType i=0; i<static_cast<IndexType>(p_target_dimension); i++)
390 embedding.col(i).array() /= sqrt(sqrt(landmarks_embedding.second(i)));
391 return TapkeeOutput(embedding,unimplementedProjectingFunction());
396 Neighbors neighbors = findNeighborsWith(kernel_distance);
401 features,current_dimension);
408 return TapkeeOutput(
project(projection_result.first,mean_vector,begin,end,features,current_dimension),projecting_function);
413 Neighbors neighbors = findNeighborsWith(kernel_distance);
418 unimplementedProjectingFunction());
423 Neighbors neighbors = findNeighborsWith(plain_distance);
428 unimplementedProjectingFunction());
433 Neighbors neighbors = findNeighborsWith(plain_distance);
438 features,current_dimension);
441 SmallestEigenvalues,eigenproblem_matrices.first,eigenproblem_matrices.second,p_target_dimension);
445 return TapkeeOutput(
project(projection_result.first,mean_vector,begin,end,features,current_dimension), projecting_function);
458 return TapkeeOutput(
project(projection_result.first,mean_vector,begin,end,features,current_dimension), projecting_function);
469 return TapkeeOutput(
project(projection_matrix,mean_vector,begin,end,features,current_dimension), projecting_function);
478 for (
IndexType i=0; i<static_cast<IndexType>(p_target_dimension); i++)
479 embedding.first.col(i).array() *= sqrt(embedding.second(i));
480 return TapkeeOutput(embedding.first, unimplementedProjectingFunction());
485 Neighbors neighbors = findNeighborsWith(kernel_distance);
490 features,current_dimension);
493 eig_matrices.first,eig_matrices.second,p_target_dimension);
497 return TapkeeOutput(
project(projection_result.first,mean_vector,begin,end,features,current_dimension),
498 projecting_function);
504 if (p_global_strategy.
is(
false))
506 neighbors = findNeighborsWith(plain_distance);
510 p_target_dimension,p_global_strategy,p_tolerance,p_n_updates,p_max_iteration), unimplementedProjectingFunction());
523 return TapkeeOutput(
project(begin,end,features,current_dimension,p_max_iteration,p_epsilon,
524 p_target_dimension, mean_vector), unimplementedProjectingFunction());
534 DenseMatrix embedding(static_cast<IndexType>(p_target_dimension),n_vectors);
536 tsne.
run(data,data.cols(),data.rows(),embedding.data(),p_target_dimension,p_perplexity,p_theta);
538 return TapkeeOutput(embedding.transpose(), unimplementedProjectingFunction());
548 Neighbors neighbors = findNeighborsWith(plain_distance);
550 manifold_sculpting_embed(begin, end, embedding, p_target_dimension, neighbors, distance, p_max_iteration, p_squishing_rate);
557 template <
class RandomAccessIterator,
class KernelCallback,
558 class DistanceCallback,
class FeaturesCallback>
560 initialize(RandomAccessIterator begin, RandomAccessIterator end,
561 KernelCallback kernel, DistanceCallback
distance, FeaturesCallback features,
565 begin,end,kernel,
distance,features,pmap,ctx);
ImplementationBase(RandomAccessIterator b, RandomAccessIterator e, KernelCallback k, DistanceCallback d, FeaturesCallback f, ParametersSet &pmap, const Context &ctx)
static const EigenMethod Dense("Dense")
Eigen library dense method (could be useful for debugging). Computes all eigenvectors thus can be ver...
static const EigendecompositionStrategy SmallestEigenvalues("Smallest eigenvalues", 1)
std::string failureMessage(const stichwort::Parameter &p) const
ScalarType distance(Callback &cb, const CoverTreePoint< RandomAccessIterator > &l, const CoverTreePoint< RandomAccessIterator > &r, ScalarType upper_bound)
TapkeeOutput embedDiffusionMap()
static const EigendecompositionStrategy SquaredLargestEigenvalues("Largest eigenvalues of squared matrix", 0)
CheckedParameter checked()
void centerMatrix(DenseMatrix &matrix)
Neighbors find_neighbors(NeighborsMethod method, const RandomAccessIterator &begin, const RandomAccessIterator &end, const Callback &callback, IndexType k, bool check_connectivity)
TapkeeOutput embedLaplacianEigenmaps()
bool operator()(T v) const
Eigen::Matrix< tapkee::ScalarType, Eigen::Dynamic, Eigen::Dynamic > DenseMatrix
dense matrix type (non-overridable)
EigendecompositionResult eigendecomposition(const EigenMethod &method, const ComputationStrategy &strategy, const EigendecompositionStrategy &eigen_strategy, const MatrixType &m, IndexType target_dimension)
Multiple implementation handler method for various eigendecomposition methods.
DenseSymmetricMatrix compute_centered_kernel_matrix(RandomAccessIterator begin, RandomAccessIterator end, KernelCallback callback)
void manifold_sculpting_embed(RandomAccessIterator begin, RandomAccessIterator end, DenseMatrix &data, IndexType target_dimension, const Neighbors &neighbors, DistanceCallback callback, IndexType max_iteration, ScalarType squishing_rate)
Parameter & satisfies(const F< Q > &cond) const
DenseSymmetricMatrix compute_covariance_matrix(RandomAccessIterator begin, RandomAccessIterator end, const DenseVector &mean, FeatureVectorCallback callback, IndexType dimension)
Namespace containing implementation of t-SNE algorithm.
TapkeeOutput embedKernelLocallyLinearEmbedding()
DenseSymmetricMatrixPair construct_lltsa_eigenproblem(SparseWeightMatrix W, RandomAccessIterator begin, RandomAccessIterator end, FeatureVectorCallback feature_vector_callback, IndexType dimension)
TapkeeOutput embedMultidimensionalScaling()
SparseWeightMatrix hessian_weight_matrix(RandomAccessIterator begin, RandomAccessIterator end, const Neighbors &neighbors, PairwiseCallback callback, const IndexType target_dimension)
DenseSymmetricMatrix compute_shortest_distances_matrix(RandomAccessIterator begin, RandomAccessIterator end, Neighbors &neighbors, DistanceCallback callback)
Computes shortest distances (so-called geodesic distances) using Dijkstra algorithm.
static Parameter create(const std::string &name, const T &value)
TAPKEE_INTERNAL_PAIR< tapkee::SparseWeightMatrix, tapkee::DenseDiagonalMatrix > Laplacian
bool operator()(T v) const
TapkeeOutput embedNeighborhoodPreservingEmbedding()
IndexType current_dimension
TapkeeOutput embedLocalityPreservingProjections()
TapkeeOutput embedHessianLocallyLinearEmbedding()
TapkeeOutput embedEmpty()
TapkeeOutput embedManifoldSculpting()
TapkeeOutput embedFactorAnalysis()
double ScalarType
default scalar value (can be overrided with TAPKEE_CUSTOM_INTERNAL_NUMTYPE define) ...
EigendecompositionResult generalized_eigendecomposition(const EigenMethod &method, const ComputationStrategy &strategy, const EigendecompositionStrategy &eigen_strategy, const LMatrixType &lhs, const RMatrixType &rhs, IndexType target_dimension)
ParameterName name() const
Parameter p_global_strategy
Parameter p_target_dimension
An exception type that is thrown when no data is given.
TAPKEE_INTERNAL_PAIR< tapkee::DenseSymmetricMatrix, tapkee::DenseSymmetricMatrix > DenseSymmetricMatrixPair
void run(tapkee::DenseMatrix &X, int N, int D, ScalarType *Y, int no_dims, ScalarType perplexity, ScalarType theta)
DenseSymmetricMatrix compute_distance_matrix(RandomAccessIterator begin, RandomAccessIterator, Landmarks &landmarks, PairwiseCallback callback)
TapkeeOutput embedPassThru()
DenseSymmetricMatrixPair construct_neighborhood_preserving_eigenproblem(SparseWeightMatrix W, RandomAccessIterator begin, RandomAccessIterator end, FeatureVectorCallback feature_vector_callback, IndexType dimension)
SparseWeightMatrix tangent_weight_matrix(RandomAccessIterator begin, RandomAccessIterator end, const Neighbors &neighbors, PairwiseCallback callback, const IndexType target_dimension, const ScalarType shift)
An exception type that is thrown when computations were cancelled.
Eigen::SparseMatrix< tapkee::ScalarType > SparseWeightMatrix
sparse weight matrix type (non-overridable)
Parameter p_max_iteration
bool operator()(T v) const
TapkeeOutput embedtDistributedStochasticNeighborEmbedding()
Parameter p_neighbors_method
DenseMatrix triangulate(RandomAccessIterator begin, RandomAccessIterator end, PairwiseCallback distance_callback, Landmarks &landmarks, DenseVector &landmark_distances_squared, EigendecompositionResult &landmarks_embedding, IndexType target_dimension)
Eigen::Matrix< tapkee::ScalarType, Eigen::Dynamic, 1 > DenseVector
dense vector type (non-overridable)
DenseMatrix dense_matrix_from_features(FeaturesCallback features, IndexType dimension, RandomAccessIterator begin, RandomAccessIterator end)
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
TapkeeOutput embedUsing(DimensionReductionMethod method)
PlainDistance< RandomAccessIterator, DistanceCallback > plain_distance
Parameter p_computation_strategy
DimensionReductionMethod
Dimension reduction methods.
TapkeeOutput embedLinearLocalTangentSpaceAlignment()
std::string failureMessage(const stichwort::Parameter &p) const
DistanceCallback distance
TapkeeOutput embedIsomap()
RandomAccessIterator begin
DenseVector compute_mean(RandomAccessIterator begin, RandomAccessIterator end, FeatureVectorCallback callback, IndexType dimension)
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)
DenseSymmetricMatrix compute_diffusion_matrix(RandomAccessIterator begin, RandomAccessIterator end, DistanceCallback callback, const ScalarType width)
Computes diffusion process matrix. Uses the following algorithm:
static const EigendecompositionStrategy LargestEigenvalues("Largest eigenvalues", 0)
A pimpl wrapper for projecting function.
bool operator()(T v) const
static tapkee::ProjectingFunction unimplementedProjectingFunction()
Return result of the library - a pair of DenseMatrix (embedding) and ProjectingFunction.
TAPKEE_INTERNAL_VECTOR< tapkee::tapkee_internal::LocalNeighbors > Neighbors
ImplementationBase< RandomAccessIterator, KernelCallback, DistanceCallback, FeaturesCallback > initialize(RandomAccessIterator begin, RandomAccessIterator end, KernelCallback kernel, DistanceCallback distance, FeaturesCallback features, stichwort::ParametersSet &pmap, const Context &ctx)
Neighbors findNeighborsWith(Distance d)
SparseWeightMatrix linear_weight_matrix(const RandomAccessIterator &begin, const RandomAccessIterator &end, const Neighbors &neighbors, PairwiseCallback callback, const ScalarType shift, const ScalarType trace_shift)
TapkeeOutput embedKernelLocalTangentSpaceAlignment()
#define tapkee_method_handle(X)
TapkeeOutput embedKernelPCA()
std::string failureMessage(const stichwort::Parameter &p) const
An exception type that is thrown when unsupported method is called.
Parameter p_check_connectivity
TAPKEE_INTERNAL_PAIR< tapkee::DenseMatrix, tapkee::DenseVector > EigendecompositionResult
TapkeeOutput embedRandomProjection()
TapkeeOutput embedLandmarkMultidimensionalScaling()
TapkeeOutput embedStochasticProximityEmbedding()
Landmarks select_landmarks_random(RandomAccessIterator begin, RandomAccessIterator end, ScalarType ratio)
std::string failureMessage(const stichwort::Parameter &p) const
TAPKEE_INTERNAL_VECTOR< tapkee::IndexType > Landmarks
KernelDistance< RandomAccessIterator, KernelCallback > kernel_distance
tapkee::DenseMatrix DenseSymmetricMatrix
dense symmetric matrix (non-overridable, currently just dense matrix, can be improved later) ...
DenseSymmetricMatrixPair construct_locality_preserving_eigenproblem(SparseWeightMatrix &L, DenseDiagonalMatrix &D, RandomAccessIterator begin, RandomAccessIterator end, FeatureVectorCallback feature_vector_callback, IndexType dimension)
Parameter p_squishing_rate
DenseMatrix gaussian_projection_matrix(IndexType target_dimension, IndexType current_dimension)
FeaturesCallback features
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)
Basic ProjectionImplementation that subtracts mean from the vector and multiplies projecting matrix w...
TapkeeOutput embedLandmarkIsomap()