6 #ifndef TAPKEE_MANIFOLD_SCULPTING_H_
7 #define TAPKEE_MANIFOLD_SCULPTING_H_
22 namespace tapkee_internal
27 const ScalarType max_number_of_iterations_without_improvement = 20;
29 const ScalarType weight_for_adjusted_point = 10.0;
30 const ScalarType learning_rate_grow_factor = 1.1;
31 const ScalarType learning_rate_shrink_factor = 0.9;
66 template<
class RandomAccessIterator,
class DistanceCallback>
68 const Neighbors& neighbors, DistanceCallback callback,
74 throw std::runtime_error(
"Wrong size");
76 sparse_triplets.reserve(k*n);
85 current_distance = callback.distance(begin[i], begin[current_neighbors[j]]);
86 average_distance += current_distance;
87 SparseTriplet triplet(i, current_neighbors[j], current_distance);
88 sparse_triplets.push_back(triplet);
91 average_distance /= (k*n);
102 sparse_triplets.reserve(k * n_vectors);
104 Neighbors most_collinear_neighbors_of_neighbors;
105 most_collinear_neighbors_of_neighbors.reserve(n_vectors);
107 for (
IndexType i = 0; i < n_vectors; ++i)
111 most_collinear_current_neighbors.reserve(k);
115 const LocalNeighbors& neighbors_of_neighbor = neighbors[current_neighbors[j]];
117 ScalarType min_cos_value = 1.0, current_cos_value;
119 most_collinear_current_neighbors.push_back(0);
123 DenseVector neighbor_to_point = data.col(i) - data.col(current_neighbors[j]);
124 DenseVector neighbor_to_its_neighbor = data.col(neighbors_of_neighbor[l])
125 - data.col(current_neighbors[j]);
126 current_cos_value = neighbor_to_point.dot(neighbor_to_its_neighbor) /
127 (neighbor_to_point.norm() *
128 neighbor_to_its_neighbor.norm());
129 if (current_cos_value < min_cos_value)
131 most_collinear_current_neighbors[j] = neighbors_of_neighbor[l];
132 min_cos_value = current_cos_value;
136 SparseTriplet triplet(i, most_collinear_current_neighbors[j], min_cos_value);
137 sparse_triplets.push_back(triplet);
139 most_collinear_neighbors_of_neighbors.push_back(most_collinear_current_neighbors);
143 most_collinear_neighbors_of_neighbors);
151 for (
IndexType i = 0; i < data.cols(); ++i)
155 average_distance += (data.col(i) - data.col(neighbors[i][j])).norm();
158 return average_distance / (k * data.cols());
172 DenseVector neighbor_to_point = data.col(index) - data.col(neighbor);
173 DenseVector neighbor_to_its_neighbor = data.col(neighbor_of_neighbor)
174 - data.col(neighbor);
175 ScalarType current_cos_value = neighbor_to_point.dot(neighbor_to_its_neighbor) /
176 (neighbor_to_point.norm() *
177 neighbor_to_its_neighbor.norm());
179 ScalarType current_distance = (data.col(index) - data.col(neighbor)).norm();
182 current_cos_value - error_func_data.
angles_matrix.coeff(index, neighbor_of_neighbor);
186 current_distance - error_func_data.
distance_matrix.coeff(index, neighbor);
192 (error_func_data.
adjusted_points.count(neighbor) == 0) ? 1 : weight_for_adjusted_point;
193 error_value += weight * (diff_cos * diff_cos + diff_distance * diff_distance);
230 for (
IndexType i = 0; i < target_dimension; ++i)
233 data(i, index) += learning_rate;
235 if (new_error >= old_error)
238 data(i, index) -= 2 * learning_rate;
242 if (new_error >= old_error)
244 data(i, index) += learning_rate;
248 old_error = new_error;
258 template <
class RandomAccessIterator,
class DistanceCallback>
261 const Neighbors& neighbors, DistanceCallback callback,
277 ScalarType no_improvement_counter = 0, normal_counter = 0;
278 ScalarType current_multiplier = squishing_rate;
279 ScalarType learning_rate = initial_average_distance;
280 ScalarType best_error = DBL_MAX, current_error, point_error;
284 while (((no_improvement_counter++ < max_number_of_iterations_without_improvement)
285 || (current_multiplier > multiplier_treshold))
286 && (normal_counter++ < max_iteration))
291 data.bottomRows(data.rows() - target_dimension) *= squishing_rate;
294 data.topRows(target_dimension) /= squishing_rate;
296 current_multiplier *= squishing_rate;
303 IndexType start_point_index = std::rand() % data.cols();
304 std::deque<IndexType> points_to_adjust;
305 points_to_adjust.push_back(start_point_index);
308 std::set<IndexType> adjusted_points;
310 while (!points_to_adjust.empty())
312 IndexType current_point_index = points_to_adjust.front();
313 points_to_adjust.pop_front();
314 if (adjusted_points.count(current_point_index) == 0)
317 distances_to_neighbors,
318 angles_and_neighbors.first,
320 angles_and_neighbors.second,
322 initial_average_distance
325 learning_rate, error_func_data, point_error);
326 current_error += point_error;
328 std::copy(neighbors[current_point_index].begin(),
329 neighbors[current_point_index].end(),
330 std::back_inserter(points_to_adjust));
331 adjusted_points.insert(current_point_index);
335 if (steps_made > data.cols())
336 learning_rate *= learning_rate_grow_factor;
338 learning_rate *= learning_rate_shrink_factor;
340 if (current_error < best_error)
342 best_error = current_error;
343 no_improvement_counter = 0;
347 data.conservativeResize(target_dimension, Eigen::NoChange);
348 data.transposeInPlace();
const Neighbors & distance_neighbors
SparseMatrixNeighborsPair angles_matrix_and_neighbors(const Neighbors &neighbors, const DenseMatrix &data)
const Neighbors & angle_neighbors
Eigen::Matrix< tapkee::ScalarType, Eigen::Dynamic, Eigen::Dynamic > DenseMatrix
dense matrix type (non-overridable)
void manifold_sculpting_embed(RandomAccessIterator begin, RandomAccessIterator end, DenseMatrix &data, IndexType target_dimension, const Neighbors &neighbors, DistanceCallback callback, IndexType max_iteration, ScalarType squishing_rate)
TAPKEE_INTERNAL_VECTOR< tapkee::tapkee_internal::SparseTriplet > SparseTriplets
const SparseMatrix & angles_matrix
double ScalarType
default scalar value (can be overrided with TAPKEE_CUSTOM_INTERNAL_NUMTYPE define) ...
ScalarType compute_error_for_point(const IndexType index, const DenseMatrix &data, const DataForErrorFunc &error_func_data)
const std::set< IndexType > & adjusted_points
Eigen::Matrix< tapkee::ScalarType, Eigen::Dynamic, 1 > DenseVector
dense vector type (non-overridable)
TAPKEE_INTERNAL_VECTOR< tapkee::IndexType > LocalNeighbors
IndexType adjust_point_at_index(const IndexType index, DenseMatrix &data, const IndexType target_dimension, const ScalarType learning_rate, const DataForErrorFunc &error_func_data, ScalarType &point_error)
SparseMatrix sparse_matrix_from_triplets(const SparseTriplets &sparse_triplets, IndexType m, IndexType n)
int IndexType
indexing type (non-overridable) set to int for compatibility with OpenMP 2.0
Eigen::Triplet< tapkee::ScalarType > SparseTriplet
Eigen::SparseMatrix< tapkee::ScalarType > SparseMatrix
sparse matrix type (non-overridable)
TAPKEE_INTERNAL_VECTOR< tapkee::tapkee_internal::LocalNeighbors > Neighbors
const SparseMatrix & distance_matrix
TAPKEE_INTERNAL_PAIR< tapkee::SparseMatrix, tapkee::tapkee_internal::Neighbors > SparseMatrixNeighborsPair
Data needed to compute error function.
const ScalarType average_distance
ScalarType average_neighbor_distance(const DenseMatrix &data, const Neighbors &neighbors)
SparseMatrix neighbors_distances_matrix(RandomAccessIterator begin, RandomAccessIterator end, const Neighbors &neighbors, DistanceCallback callback, ScalarType &average_distance)