Tapkee
methods.hpp
Go to the documentation of this file.
00001 /* This software is distributed under BSD 3-clause license (see LICENSE file).
00002  *
00003  * Copyright (c) 2012-2013 Sergey Lisitsyn, Fernando Iglesias
00004  */
00005 
00006 #ifndef TAPKEE_METHODS_H_
00007 #define TAPKEE_METHODS_H_
00008 
00009 /* Tapkee includes */
00010 #include <tapkee/defines.hpp>
00011 #include <tapkee/utils/naming.hpp>
00012 #include <tapkee/utils/time.hpp>
00013 #include <tapkee/utils/logging.hpp>
00014 #include <tapkee/utils/conditional_select.hpp>
00015 #include <tapkee/utils/features.hpp>
00016 #include <tapkee/parameters/defaults.hpp>
00017 #include <tapkee/parameters/context.hpp>
00018 #include <tapkee/routines/locally_linear.hpp>
00019 #include <tapkee/routines/eigendecomposition.hpp>
00020 #include <tapkee/routines/generalized_eigendecomposition.hpp>
00021 #include <tapkee/routines/multidimensional_scaling.hpp>
00022 #include <tapkee/routines/diffusion_maps.hpp>
00023 #include <tapkee/routines/laplacian_eigenmaps.hpp>
00024 #include <tapkee/routines/isomap.hpp>
00025 #include <tapkee/routines/pca.hpp>
00026 #include <tapkee/routines/random_projection.hpp>
00027 #include <tapkee/routines/spe.hpp>
00028 #include <tapkee/routines/fa.hpp>
00029 #include <tapkee/routines/manifold_sculpting.hpp>
00030 #include <tapkee/neighbors/neighbors.hpp>
00031 #include <tapkee/external/barnes_hut_sne/tsne.hpp>
00032 /* End of Tapkee includes */
00033 
00034 namespace tapkee
00035 {
00037 namespace tapkee_internal
00038 {
00039 
00040 template <class RandomAccessIterator, class KernelCallback,
00041           class DistanceCallback, class FeaturesCallback>
00042 class ImplementationBase
00043 {
00044 public:
00045 
00046     ImplementationBase(RandomAccessIterator b, RandomAccessIterator e,
00047                        KernelCallback k, DistanceCallback d, FeaturesCallback f,
00048                        ParametersSet& pmap, const Context& ctx) : 
00049         parameters(pmap), context(ctx), kernel(k), distance(d), features(f),
00050         plain_distance(PlainDistance<RandomAccessIterator,DistanceCallback>(distance)),
00051         kernel_distance(KernelDistance<RandomAccessIterator,KernelCallback>(kernel)),
00052         begin(b), end(e), 
00053         eigen_method(), neighbors_method(), eigenshift(), traceshift(),
00054         check_connectivity(), n_neighbors(), width(), timesteps(),
00055         ratio(), max_iteration(), tolerance(), n_updates(), perplexity(), 
00056         theta(), squishing_rate(), global_strategy(), epsilon(), target_dimension(),
00057         n_vectors(0), current_dimension(0)
00058     {
00059         n_vectors = (end-begin);
00060 
00061         target_dimension = parameters(keywords::target_dimension);
00062         n_neighbors = parameters(keywords::num_neighbors).checked().positive();
00063         
00064         if (n_vectors > 0)
00065         {
00066             target_dimension.checked()
00067                 .inRange(static_cast<IndexType>(1),static_cast<IndexType>(n_vectors));
00068             n_neighbors.checked()
00069                 .inRange(static_cast<IndexType>(3),static_cast<IndexType>(n_vectors));
00070         }
00071 
00072         eigen_method = parameters(keywords::eigen_method);
00073         neighbors_method = parameters(keywords::neighbors_method);
00074         check_connectivity = parameters(keywords::check_connectivity);
00075         width = parameters(keywords::gaussian_kernel_width).checked().positive();
00076         timesteps = parameters(keywords::diffusion_map_timesteps).checked().positive();
00077         eigenshift = parameters(keywords::nullspace_shift);
00078         traceshift = parameters(keywords::klle_shift);
00079         max_iteration = parameters(keywords::max_iteration);
00080         tolerance = parameters(keywords::spe_tolerance).checked().positive();
00081         n_updates = parameters(keywords::spe_num_updates).checked().positive();
00082         theta = parameters(keywords::sne_theta).checked().nonNegative();
00083         squishing_rate = parameters(keywords::squishing_rate);
00084         global_strategy = parameters(keywords::spe_global_strategy);
00085         epsilon = parameters(keywords::fa_epsilon).checked().nonNegative();
00086         perplexity = parameters(keywords::sne_perplexity).checked().nonNegative();
00087         ratio = parameters(keywords::landmark_ratio);
00088 
00089         if (!is_dummy<FeaturesCallback>::value)
00090         {
00091             current_dimension = features.dimension();
00092         }
00093         else
00094         {
00095             current_dimension = 0;
00096         }
00097     }
00098 
00099     TapkeeOutput embedUsing(DimensionReductionMethod method)
00100     {
00101         if (context.is_cancelled()) 
00102             throw cancelled_exception();
00103 
00104         using std::mem_fun_ref_t;
00105         using std::mem_fun_ref;
00106         typedef std::mem_fun_ref_t<TapkeeOutput,ImplementationBase> ImplRef;
00107 
00108 #define tapkee_method_handle(X)                                                                 \
00109         case X:                                                                                 \
00110         {                                                                                       \
00111             timed_context tctx__("[+] embedding with " # X);                                    \
00112             ImplRef ref = conditional_select<                                                   \
00113                 ((!MethodTraits<X>::needs_kernel)   || (!is_dummy<KernelCallback>::value))   && \
00114                 ((!MethodTraits<X>::needs_distance) || (!is_dummy<DistanceCallback>::value)) && \
00115                 ((!MethodTraits<X>::needs_features) || (!is_dummy<FeaturesCallback>::value)),   \
00116                     ImplRef>()(mem_fun_ref(&ImplementationBase::embed##X),                      \
00117                                mem_fun_ref(&ImplementationBase::embedEmpty));                   \
00118             return ref(*this);                                                                  \
00119         }                                                                                       \
00120         break                                                                                   \
00121 
00122         switch (method)
00123         {
00124             tapkee_method_handle(KernelLocallyLinearEmbedding);
00125             tapkee_method_handle(KernelLocalTangentSpaceAlignment);
00126             tapkee_method_handle(DiffusionMap);
00127             tapkee_method_handle(MultidimensionalScaling);
00128             tapkee_method_handle(LandmarkMultidimensionalScaling);
00129             tapkee_method_handle(Isomap);
00130             tapkee_method_handle(LandmarkIsomap);
00131             tapkee_method_handle(NeighborhoodPreservingEmbedding);
00132             tapkee_method_handle(LinearLocalTangentSpaceAlignment);
00133             tapkee_method_handle(HessianLocallyLinearEmbedding);
00134             tapkee_method_handle(LaplacianEigenmaps);
00135             tapkee_method_handle(LocalityPreservingProjections);
00136             tapkee_method_handle(PCA);
00137             tapkee_method_handle(KernelPCA);
00138             tapkee_method_handle(RandomProjection);
00139             tapkee_method_handle(StochasticProximityEmbedding);
00140             tapkee_method_handle(PassThru);
00141             tapkee_method_handle(FactorAnalysis);
00142             tapkee_method_handle(tDistributedStochasticNeighborEmbedding);
00143             tapkee_method_handle(ManifoldSculpting);
00144         }
00145 #undef tapkee_method_handle
00146         return TapkeeOutput();
00147     }
00148 
00149 private:
00150 
00151     static const IndexType SkipOneEigenvalue = 1;
00152     static const IndexType SkipNoEigenvalues = 0;
00153 
00154     ParametersSet parameters;
00155     Context context;
00156     KernelCallback kernel;
00157     DistanceCallback distance;
00158     FeaturesCallback features;
00159     PlainDistance<RandomAccessIterator,DistanceCallback> plain_distance;
00160     KernelDistance<RandomAccessIterator,KernelCallback> kernel_distance;
00161 
00162     RandomAccessIterator begin;
00163     RandomAccessIterator end;
00164 
00165     Parameter eigen_method;
00166     Parameter neighbors_method;
00167     Parameter eigenshift;
00168     Parameter traceshift;
00169     Parameter check_connectivity;
00170     Parameter n_neighbors;
00171     Parameter width;
00172     Parameter timesteps;
00173     Parameter ratio;
00174     Parameter max_iteration;
00175     Parameter tolerance;
00176     Parameter n_updates;
00177     Parameter perplexity;
00178     Parameter theta;
00179     Parameter squishing_rate;   
00180     Parameter global_strategy;
00181     Parameter epsilon;
00182     Parameter target_dimension;
00183 
00184     IndexType n_vectors;
00185     IndexType current_dimension;
00186 
00187     template<class Distance>
00188     Neighbors findNeighborsWith(Distance d)
00189     {
00190         return find_neighbors(neighbors_method,begin,end,d,n_neighbors,check_connectivity);
00191     }
00192 
00193     static tapkee::ProjectingFunction unimplementedProjectingFunction() 
00194     {
00195         return tapkee::ProjectingFunction();
00196     }
00197 
00198     TapkeeOutput embedEmpty()
00199     {
00200         throw unsupported_method_error("Some callback is missed");
00201         return TapkeeOutput();
00202     }
00203 
00204     TapkeeOutput embedKernelLocallyLinearEmbedding()
00205     {
00206         Neighbors neighbors = findNeighborsWith(kernel_distance);
00207         SparseWeightMatrix weight_matrix =
00208             linear_weight_matrix(begin,end,neighbors,kernel,eigenshift,traceshift);
00209         DenseMatrix embedding = 
00210             eigendecomposition<SparseWeightMatrix,SparseInverseMatrixOperation>(eigen_method,
00211                 weight_matrix,target_dimension,SkipOneEigenvalue).first;
00212 
00213         return TapkeeOutput(embedding, unimplementedProjectingFunction());
00214     }
00215 
00216     TapkeeOutput embedKernelLocalTangentSpaceAlignment()
00217     {
00218         Neighbors neighbors = findNeighborsWith(kernel_distance);
00219         SparseWeightMatrix weight_matrix = 
00220             tangent_weight_matrix(begin,end,neighbors,kernel,target_dimension,eigenshift);
00221         DenseMatrix embedding =
00222             eigendecomposition<SparseWeightMatrix,SparseInverseMatrixOperation>(eigen_method,
00223                 weight_matrix,target_dimension,SkipOneEigenvalue).first;
00224 
00225         return TapkeeOutput(embedding, unimplementedProjectingFunction());
00226     }
00227 
00228     TapkeeOutput embedDiffusionMap()
00229     {
00230         #ifdef TAPKEE_GPU
00231             #define DM_MATRIX_OP GPUDenseImplicitSquareMatrixOperation
00232         #else
00233             #define DM_MATRIX_OP DenseImplicitSquareSymmetricMatrixOperation
00234         #endif
00235 
00236         DenseSymmetricMatrix diffusion_matrix =
00237             compute_diffusion_matrix(begin,end,distance,timesteps,width);
00238         DenseMatrix embedding =
00239             eigendecomposition<DenseSymmetricMatrix,DM_MATRIX_OP>(eigen_method,diffusion_matrix,
00240                 target_dimension,SkipNoEigenvalues).first;
00241 
00242         return TapkeeOutput(embedding, unimplementedProjectingFunction());
00243 
00244         #undef DM_MATRIX_OP
00245     }
00246 
00247     TapkeeOutput embedMultidimensionalScaling()
00248     {
00249         #ifdef TAPKEE_GPU
00250             #define MDS_MATRIX_OP GPUDenseImplicitSquareMatrixOperation
00251         #else
00252             #define MDS_MATRIX_OP DenseMatrixOperation
00253         #endif
00254 
00255         DenseSymmetricMatrix distance_matrix = compute_distance_matrix(begin,end,distance);
00256         centerMatrix(distance_matrix);
00257         distance_matrix.array() *= -0.5;
00258         EigendecompositionResult embedding = 
00259             eigendecomposition<DenseSymmetricMatrix,MDS_MATRIX_OP>(eigen_method,
00260                 distance_matrix,target_dimension,SkipNoEigenvalues);
00261 
00262         for (IndexType i=0; i<static_cast<IndexType>(target_dimension); i++)
00263             embedding.first.col(i).array() *= sqrt(embedding.second(i));
00264         return TapkeeOutput(embedding.first, unimplementedProjectingFunction());
00265         #undef MDS_MATRIX_OP
00266     }
00267 
00268     TapkeeOutput embedLandmarkMultidimensionalScaling()
00269     {
00270         ratio.checked()
00271             .inClosedRange(static_cast<ScalarType>(3.0/n_vectors),
00272                            static_cast<ScalarType>(1.0));
00273 
00274         Landmarks landmarks = 
00275             select_landmarks_random(begin,end,ratio);
00276         DenseSymmetricMatrix distance_matrix = 
00277             compute_distance_matrix(begin,end,landmarks,distance);
00278         DenseVector landmark_distances_squared = distance_matrix.colwise().mean();
00279         centerMatrix(distance_matrix);
00280         distance_matrix.array() *= -0.5;
00281         EigendecompositionResult landmarks_embedding = 
00282             eigendecomposition<DenseSymmetricMatrix,DenseMatrixOperation>(eigen_method,
00283                 distance_matrix,target_dimension,SkipNoEigenvalues);
00284         for (IndexType i=0; i<static_cast<IndexType>(target_dimension); i++)
00285             landmarks_embedding.first.col(i).array() *= sqrt(landmarks_embedding.second(i));
00286         return TapkeeOutput(triangulate(begin,end,distance,landmarks,
00287             landmark_distances_squared,landmarks_embedding,target_dimension), unimplementedProjectingFunction());
00288     }
00289 
00290     TapkeeOutput embedIsomap()
00291     {
00292         Neighbors neighbors = findNeighborsWith(plain_distance);
00293         DenseSymmetricMatrix shortest_distances_matrix = 
00294             compute_shortest_distances_matrix(begin,end,neighbors,distance);
00295         shortest_distances_matrix = shortest_distances_matrix.array().square();
00296         centerMatrix(shortest_distances_matrix);
00297         shortest_distances_matrix.array() *= -0.5;
00298         
00299         EigendecompositionResult embedding = 
00300             eigendecomposition<DenseSymmetricMatrix,DenseMatrixOperation>(eigen_method,
00301                 shortest_distances_matrix,target_dimension,SkipNoEigenvalues);
00302 
00303         for (IndexType i=0; i<static_cast<IndexType>(target_dimension); i++)
00304             embedding.first.col(i).array() *= sqrt(embedding.second(i));
00305 
00306         return TapkeeOutput(embedding.first, unimplementedProjectingFunction());
00307     }
00308 
00309     TapkeeOutput embedLandmarkIsomap()
00310     {
00311         ratio.checked()
00312             .inClosedRange(static_cast<ScalarType>(3.0/n_vectors),
00313                            static_cast<ScalarType>(1.0));
00314 
00315         Neighbors neighbors = findNeighborsWith(plain_distance);
00316         Landmarks landmarks = 
00317             select_landmarks_random(begin,end,ratio);
00318         DenseMatrix distance_matrix = 
00319             compute_shortest_distances_matrix(begin,end,landmarks,neighbors,distance);
00320         distance_matrix = distance_matrix.array().square();
00321         
00322         DenseVector col_means = distance_matrix.colwise().mean();
00323         DenseVector row_means = distance_matrix.rowwise().mean();
00324         ScalarType grand_mean = distance_matrix.mean();
00325         distance_matrix.array() += grand_mean;
00326         distance_matrix.colwise() -= row_means;
00327         distance_matrix.rowwise() -= col_means.transpose();
00328         distance_matrix.array() *= -0.5;
00329 
00330         EigendecompositionResult landmarks_embedding;
00331         
00332         if (eigen_method.is(Dense))
00333         {
00334             DenseMatrix distance_matrix_sym = distance_matrix*distance_matrix.transpose();
00335             landmarks_embedding = eigendecomposition<DenseSymmetricMatrix,DenseImplicitSquareMatrixOperation>
00336                 (eigen_method,distance_matrix_sym,target_dimension,SkipNoEigenvalues);
00337         }
00338         else 
00339         {
00340             landmarks_embedding = eigendecomposition<DenseSymmetricMatrix,DenseImplicitSquareMatrixOperation>
00341                 (eigen_method,distance_matrix,target_dimension,SkipNoEigenvalues);
00342         }
00343 
00344         DenseMatrix embedding = distance_matrix.transpose()*landmarks_embedding.first;
00345 
00346         for (IndexType i=0; i<static_cast<IndexType>(target_dimension); i++)
00347             embedding.col(i).array() /= sqrt(sqrt(landmarks_embedding.second(i)));
00348         return TapkeeOutput(embedding,unimplementedProjectingFunction());
00349     }
00350 
00351     TapkeeOutput embedNeighborhoodPreservingEmbedding()
00352     {
00353         Neighbors neighbors = findNeighborsWith(kernel_distance);
00354         SparseWeightMatrix weight_matrix = 
00355             linear_weight_matrix(begin,end,neighbors,kernel,eigenshift,traceshift);
00356         DenseSymmetricMatrixPair eig_matrices =
00357             construct_neighborhood_preserving_eigenproblem(weight_matrix,begin,end,
00358                 features,current_dimension);
00359         EigendecompositionResult projection_result = 
00360             generalized_eigendecomposition<DenseSymmetricMatrix,DenseSymmetricMatrix,DenseInverseMatrixOperation>(
00361                 eigen_method,eig_matrices.first,eig_matrices.second,target_dimension,SkipNoEigenvalues);
00362         DenseVector mean_vector = 
00363             compute_mean(begin,end,features,current_dimension);
00364         tapkee::ProjectingFunction projecting_function(new tapkee::MatrixProjectionImplementation(projection_result.first,mean_vector));
00365         return TapkeeOutput(project(projection_result.first,mean_vector,begin,end,features,current_dimension),projecting_function);
00366     }
00367 
00368     TapkeeOutput embedHessianLocallyLinearEmbedding()
00369     {
00370         Neighbors neighbors = findNeighborsWith(kernel_distance);
00371         SparseWeightMatrix weight_matrix =
00372             hessian_weight_matrix(begin,end,neighbors,kernel,target_dimension);
00373         return TapkeeOutput(eigendecomposition<SparseWeightMatrix,SparseInverseMatrixOperation>(eigen_method,
00374             weight_matrix,target_dimension,SkipOneEigenvalue).first, unimplementedProjectingFunction());
00375     }
00376 
00377     TapkeeOutput embedLaplacianEigenmaps()
00378     {
00379         Neighbors neighbors = findNeighborsWith(plain_distance);
00380         Laplacian laplacian = 
00381             compute_laplacian(begin,end,neighbors,distance,width);
00382         return TapkeeOutput(generalized_eigendecomposition<SparseWeightMatrix,DenseDiagonalMatrix,SparseInverseMatrixOperation>(
00383             eigen_method,laplacian.first,laplacian.second,target_dimension,SkipOneEigenvalue).first, unimplementedProjectingFunction());
00384     }
00385 
00386     TapkeeOutput embedLocalityPreservingProjections()
00387     {
00388         Neighbors neighbors = findNeighborsWith(plain_distance);
00389         Laplacian laplacian = 
00390             compute_laplacian(begin,end,neighbors,distance,width);
00391         DenseSymmetricMatrixPair eigenproblem_matrices =
00392             construct_locality_preserving_eigenproblem(laplacian.first,laplacian.second,begin,end,
00393                     features,current_dimension);
00394         EigendecompositionResult projection_result = 
00395             generalized_eigendecomposition<DenseSymmetricMatrix,DenseSymmetricMatrix,DenseInverseMatrixOperation>(
00396                 eigen_method,eigenproblem_matrices.first,eigenproblem_matrices.second,target_dimension,SkipNoEigenvalues);
00397         DenseVector mean_vector = 
00398             compute_mean(begin,end,features,current_dimension);
00399         tapkee::ProjectingFunction projecting_function(new tapkee::MatrixProjectionImplementation(projection_result.first,mean_vector));
00400         return TapkeeOutput(project(projection_result.first,mean_vector,begin,end,features,current_dimension), projecting_function);
00401     }
00402 
00403     TapkeeOutput embedPCA()
00404     {
00405         DenseVector mean_vector = 
00406             compute_mean(begin,end,features,current_dimension);
00407         DenseSymmetricMatrix centered_covariance_matrix = 
00408             compute_covariance_matrix(begin,end,mean_vector,features,current_dimension);
00409         EigendecompositionResult projection_result = 
00410             eigendecomposition<DenseSymmetricMatrix,DenseMatrixOperation>(eigen_method,centered_covariance_matrix,target_dimension,SkipNoEigenvalues);
00411         tapkee::ProjectingFunction projecting_function(new tapkee::MatrixProjectionImplementation(projection_result.first,mean_vector));
00412         return TapkeeOutput(project(projection_result.first,mean_vector,begin,end,features,current_dimension), projecting_function);
00413     }
00414 
00415     TapkeeOutput embedRandomProjection()
00416     {
00417         DenseMatrix projection_matrix = 
00418             gaussian_projection_matrix(current_dimension, target_dimension);
00419         DenseVector mean_vector = 
00420             compute_mean(begin,end,features,current_dimension);
00421 
00422         tapkee::ProjectingFunction projecting_function(new tapkee::MatrixProjectionImplementation(projection_matrix,mean_vector));
00423         return TapkeeOutput(project(projection_matrix,mean_vector,begin,end,features,current_dimension), projecting_function);
00424     }
00425 
00426     TapkeeOutput embedKernelPCA()
00427     {
00428         DenseSymmetricMatrix centered_kernel_matrix = 
00429             compute_centered_kernel_matrix(begin,end,kernel);
00430         EigendecompositionResult embedding = eigendecomposition<DenseSymmetricMatrix,DenseMatrixOperation>(eigen_method,
00431             centered_kernel_matrix,target_dimension,SkipNoEigenvalues);
00432         for (IndexType i=0; i<static_cast<IndexType>(target_dimension); i++)
00433             embedding.first.col(i).array() *= sqrt(embedding.second(i));
00434         return TapkeeOutput(embedding.first, unimplementedProjectingFunction());
00435     }
00436 
00437     TapkeeOutput embedLinearLocalTangentSpaceAlignment()
00438     {
00439         Neighbors neighbors = findNeighborsWith(kernel_distance);
00440         SparseWeightMatrix weight_matrix = 
00441             tangent_weight_matrix(begin,end,neighbors,kernel,target_dimension,eigenshift);
00442         DenseSymmetricMatrixPair eig_matrices =
00443             construct_lltsa_eigenproblem(weight_matrix,begin,end,
00444                 features,current_dimension);
00445         EigendecompositionResult projection_result = 
00446             generalized_eigendecomposition<DenseSymmetricMatrix,DenseSymmetricMatrix,DenseInverseMatrixOperation>(
00447                 eigen_method,eig_matrices.first,eig_matrices.second,target_dimension,SkipNoEigenvalues);
00448         DenseVector mean_vector = 
00449             compute_mean(begin,end,features,current_dimension);
00450         tapkee::ProjectingFunction projecting_function(new tapkee::MatrixProjectionImplementation(projection_result.first,mean_vector));
00451         return TapkeeOutput(project(projection_result.first,mean_vector,begin,end,features,current_dimension),
00452                 projecting_function);
00453     }
00454 
00455     TapkeeOutput embedStochasticProximityEmbedding()
00456     {
00457         Neighbors neighbors;
00458         if (global_strategy.is(false))
00459         {
00460             neighbors = findNeighborsWith(plain_distance);
00461         }
00462 
00463         return TapkeeOutput(spe_embedding(begin,end,distance,neighbors,
00464                 target_dimension,global_strategy,tolerance,n_updates,max_iteration), unimplementedProjectingFunction());
00465     }
00466 
00467     TapkeeOutput embedPassThru()
00468     {
00469         DenseMatrix feature_matrix =
00470             dense_matrix_from_features(features, current_dimension, begin, end);
00471         return TapkeeOutput(feature_matrix.transpose(),tapkee::ProjectingFunction());
00472     }
00473 
00474     TapkeeOutput embedFactorAnalysis()
00475     {
00476         DenseVector mean_vector = compute_mean(begin,end,features,current_dimension);
00477         return TapkeeOutput(project(begin,end,features,current_dimension,max_iteration,epsilon,
00478                                     target_dimension, mean_vector), unimplementedProjectingFunction());
00479     }
00480 
00481     TapkeeOutput embedtDistributedStochasticNeighborEmbedding()
00482     {
00483         perplexity.checked()
00484             .inClosedRange(static_cast<ScalarType>(0.0),
00485                            static_cast<ScalarType>((n_vectors-1)/3.0));
00486 
00487         DenseMatrix data = 
00488             dense_matrix_from_features(features, current_dimension, begin, end);
00489 
00490         DenseMatrix embedding(static_cast<IndexType>(target_dimension),n_vectors);
00491         tsne::TSNE tsne;
00492         tsne.run(data.data(),n_vectors,current_dimension,embedding.data(),target_dimension,perplexity,theta);
00493 
00494         return TapkeeOutput(embedding.transpose(), unimplementedProjectingFunction());
00495     }
00496 
00497     TapkeeOutput embedManifoldSculpting()
00498     {
00499         squishing_rate.checked()
00500             .inRange(static_cast<ScalarType>(0.0),
00501                      static_cast<ScalarType>(1.0));
00502 
00503         DenseMatrix embedding =
00504             dense_matrix_from_features(features, current_dimension, begin, end);
00505 
00506         Neighbors neighbors = findNeighborsWith(plain_distance);
00507 
00508         manifold_sculpting_embed(begin, end, embedding, target_dimension, neighbors, distance, max_iteration, squishing_rate);
00509 
00510         return TapkeeOutput(embedding, tapkee::ProjectingFunction());   
00511     }
00512 
00513 };
00514 
00515 template <class RandomAccessIterator, class KernelCallback,
00516           class DistanceCallback, class FeaturesCallback>
00517 ImplementationBase<RandomAccessIterator,KernelCallback,DistanceCallback,FeaturesCallback>
00518     initialize(RandomAccessIterator begin, RandomAccessIterator end,
00519                KernelCallback kernel, DistanceCallback distance, FeaturesCallback features,
00520                ParametersSet& pmap, const Context& ctx)
00521 {
00522     return ImplementationBase<RandomAccessIterator,KernelCallback,DistanceCallback,FeaturesCallback>(
00523             begin,end,kernel,distance,features,pmap,ctx);
00524 }
00525 
00526 } // End of namespace tapkee_internal
00527 } // End of namespace tapkee
00528 
00529 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines