Tapkee
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
laplacian_eigenmaps.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
4  */
5 
6 #ifndef TAPKEE_LaplacianEigenmaps_H_
7 #define TAPKEE_LaplacianEigenmaps_H_
8 
9 /* Tapkee includes */
10 #include <tapkee/defines.hpp>
11 #include <tapkee/utils/time.hpp>
12 /* End of Tapkee includes */
13 
14 namespace tapkee
15 {
16 namespace tapkee_internal
17 {
18 
38 template<class RandomAccessIterator, class DistanceCallback>
39 Laplacian compute_laplacian(RandomAccessIterator begin,
40  RandomAccessIterator end,const Neighbors& neighbors,
41  DistanceCallback callback, ScalarType width)
42 {
43  SparseTriplets sparse_triplets;
44 
45  timed_context context("Laplacian computation");
46  const IndexType k = neighbors[0].size();
47  sparse_triplets.reserve((k+1)*(end-begin));
48 
49  DenseVector D = DenseVector::Zero(end-begin);
50  for (RandomAccessIterator iter=begin; iter!=end; ++iter)
51  {
52  const LocalNeighbors& current_neighbors = neighbors[iter-begin];
53 
54  for (IndexType i=0; i<k; ++i)
55  {
56  ScalarType distance = callback.distance(*iter,begin[current_neighbors[i]]);
57  ScalarType heat = exp(-distance*distance/width);
58  D(iter-begin) += heat;
59  D(current_neighbors[i]) += heat;
60  sparse_triplets.push_back(SparseTriplet(current_neighbors[i],(iter-begin),-heat));
61  sparse_triplets.push_back(SparseTriplet((iter-begin),current_neighbors[i],-heat));
62  }
63  }
64  for (IndexType i=0; i<static_cast<IndexType>(end-begin); ++i)
65  sparse_triplets.push_back(SparseTriplet(i,i,D(i)));
66 
67 #ifdef EIGEN_YES_I_KNOW_SPARSE_MODULE_IS_NOT_STABLE_YET
68  Eigen::DynamicSparseMatrix<ScalarType> dynamic_weight_matrix(end-begin,end-begin);
69  dynamic_weight_matrix.reserve(sparse_triplets.size());
70  for (SparseTriplets::const_iterator it=sparse_triplets.begin(); it!=sparse_triplets.end(); ++it)
71  dynamic_weight_matrix.coeffRef(it->col(),it->row()) += it->value();
72  SparseWeightMatrix weight_matrix(dynamic_weight_matrix);
73 #else
74  SparseWeightMatrix weight_matrix(end-begin,end-begin);
75  weight_matrix.setFromTriplets(sparse_triplets.begin(),sparse_triplets.end());
76 #endif
77 
78  return Laplacian(weight_matrix,DenseDiagonalMatrix(D));
79 }
80 
81 template<class RandomAccessIterator, class FeatureVectorCallback>
83  DenseDiagonalMatrix& D, RandomAccessIterator begin, RandomAccessIterator end, FeatureVectorCallback feature_vector_callback,
84  IndexType dimension)
85 {
86  timed_context context("Constructing LPP eigenproblem");
87 
88  DenseSymmetricMatrix lhs = DenseSymmetricMatrix::Zero(dimension,dimension);
89  DenseSymmetricMatrix rhs = DenseSymmetricMatrix::Zero(dimension,dimension);
90 
91  DenseVector rank_update_vector_i(dimension);
92  DenseVector rank_update_vector_j(dimension);
93  for (RandomAccessIterator iter=begin; iter!=end; ++iter)
94  {
95  feature_vector_callback.vector(*iter,rank_update_vector_i);
96  rhs.selfadjointView<Eigen::Upper>().rankUpdate(rank_update_vector_i,D.diagonal()(iter-begin));
97  }
98 
99  for (int i=0; i<L.outerSize(); ++i)
100  {
101  for (SparseWeightMatrix::InnerIterator it(L,i); it; ++it)
102  {
103  feature_vector_callback.vector(begin[it.row()],rank_update_vector_i);
104  feature_vector_callback.vector(begin[it.col()],rank_update_vector_j);
105  lhs.selfadjointView<Eigen::Upper>().rankUpdate(rank_update_vector_i, rank_update_vector_j, it.value());
106  }
107  }
108 
109  return DenseSymmetricMatrixPair(lhs,rhs);
110 }
111 
112 }
113 }
114 
115 #endif
Eigen::DiagonalMatrix< tapkee::ScalarType, Eigen::Dynamic > DenseDiagonalMatrix
dense diagonal matrix
Definition: types.hpp:27
ScalarType distance(Callback &cb, const CoverTreePoint< RandomAccessIterator > &l, const CoverTreePoint< RandomAccessIterator > &r, ScalarType upper_bound)
TAPKEE_INTERNAL_VECTOR< tapkee::tapkee_internal::SparseTriplet > SparseTriplets
Definition: synonyms.hpp:38
TAPKEE_INTERNAL_PAIR< tapkee::SparseWeightMatrix, tapkee::DenseDiagonalMatrix > Laplacian
Definition: synonyms.hpp:43
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
Eigen::SparseMatrix< tapkee::ScalarType > SparseWeightMatrix
sparse weight matrix type (non-overridable)
Definition: types.hpp:29
Eigen::Matrix< tapkee::ScalarType, Eigen::Dynamic, 1 > DenseVector
dense vector type (non-overridable)
Definition: types.hpp:21
TAPKEE_INTERNAL_VECTOR< tapkee::IndexType > LocalNeighbors
Definition: synonyms.hpp:39
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
Eigen::Triplet< tapkee::ScalarType > SparseTriplet
Definition: synonyms.hpp:35
TAPKEE_INTERNAL_VECTOR< tapkee::tapkee_internal::LocalNeighbors > Neighbors
Definition: synonyms.hpp:40
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)