From e9be8af2ccb55b676e4e365376ac92685ad3b11e Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Tue, 2 Apr 2024 14:51:50 -0400 Subject: [PATCH 01/11] Clean up MakeAlias() and put it in only one place. --- doc/user/core.md | 80 +++++++++++- src/mlpack/core/math/make_alias.hpp | 119 ++++++++---------- src/mlpack/core/math/math.hpp | 1 + src/mlpack/core/math/unwrap_alias.hpp | 38 ++++++ src/mlpack/methods/ann/ann.hpp | 1 - src/mlpack/methods/ann/ffn_impl.hpp | 2 - src/mlpack/methods/ann/layer/add_merge.hpp | 1 - src/mlpack/methods/ann/layer/layer.hpp | 2 +- src/mlpack/methods/ann/layer/multi_layer.hpp | 1 - src/mlpack/methods/ann/make_alias.hpp | 58 --------- .../bias_svd/bias_svd_function_impl.hpp | 5 +- .../linear_svm/linear_svm_function_impl.hpp | 5 +- .../methods/lmnn/lmnn_function_impl.hpp | 13 +- .../logistic_regression_function_impl.hpp | 14 ++- .../nca/nca_softmax_error_function_impl.hpp | 11 +- .../regularized_svd_function_impl.hpp | 5 +- .../softmax_regression_function_impl.hpp | 5 +- .../svdplusplus/svdplusplus_function_impl.hpp | 5 +- 18 files changed, 202 insertions(+), 164 deletions(-) create mode 100644 src/mlpack/core/math/unwrap_alias.hpp delete mode 100644 src/mlpack/methods/ann/make_alias.hpp diff --git a/doc/user/core.md b/doc/user/core.md index e3ec6a9122..4fd8c02da8 100644 --- a/doc/user/core.md +++ b/doc/user/core.md @@ -15,6 +15,9 @@ classes, each of which are documented on this page. mlpack provides a number of additional mathematical utility classes and functions on top of Armadillo. + * [Aliases](#aliases): utilities to create and manage aliases (`MakeAlias()`, + `ClearAlias()`, `UnwrapAlias()`). + * [`Range`](#range): simple mathematical range (i.e. `[0, 3]`) * [`ColumnCovariance()`](#columncovariance): compute covariance of @@ -31,9 +34,6 @@ functions on top of Armadillo. * [Logarithmic utilities](#logarithmic-utilities): `LogAdd()`, `AccuLog()`, `LogSumExp()`, `LogSumExpT()`. - - * [`MultiplyCube2Cube()`](#multiplycube2cube): multiply each slice in a cube by each slice in another cube * [`MultiplyMat2Cube()`](#multiplymat2cube): multiply a matrix by each slice in a cube * [`MultiplyCube2Mat()`](#multiplycube2mat): multiply each slice in a cube by a matrix @@ -47,6 +47,80 @@ functions on top of Armadillo. --- +### Aliases + +Aliases are matrix, vector, or cube objects that share memory with another +matrix, vector, or cube. They are often used internally inside of mlpack to +avoid copies. + +***Important caveats about aliases***: + + - An alias represents the same memory block as the input. As such, changes to + the alias object will also be reflected in the original object. + + - The `MakeAlias()` function is not guaranteed to return an alias; it only + returns an alias *if possible*, and makes a copy otherwise. + + - If `mat` goes out of scope or is destructed, then `a` ***becomes invalid***. + _You are responsible for ensuring an invalid alias is not used!_ + +--- + + * `MakeAlias(a, mat, rows, cols, strict=true)` + - Make `a` into an alias of `mat` with the given size. + - If `strict` is `true`, the size of `a` cannot be changed. + - `mat` and `a` should have the same matrix type (e.g. `arma::mat`, + `arma::fmat`, `arma::sp_mat`). + - If an alias cannot be created, the matrix will be copied. Sparse types + cannot have aliases and will be copied. + + * `MakeAlias(a, cube, rows, cols, slices, strict=true)` + - Make `a` into an alias of `cube` with the given size. + - If `strict` is `true`, the size of `a` cannot be changed. + - `cube` and `a` should have the same matrix type (e.g. `arma::cube`, + `arma::fcube`). + - If an alias cannot be created, the matrix will be copied. + + * `MakeAlias(a, memptr, rows, cols, strict=true)` + - Make `a` into an alias of the memory block starting at `memptr` of size + `rows` by `cols`. + - The memory at `memptr` should be arranged in a [column-major + ordering](matrices.md#representing-data-in-mlpack). + - If `strict` is `true`, the size of `a` cannot be changed. + - `a` should be a dense matrix type (e.g. `arma::mat`, `arma::fmat`), and + `memptr` should be a non-const pointer of the matrix's element type (e.g. + `double*`, `float*`). + + * `MakeAlias(a, memptr, rows, cols, slices, strict=true)` + - Make `a` into an alias of the memory block starting at `memptr` of size + `rows` by `cols` by `slices`. + - The memory at `memptr` should be arranged in a [column-major + ordering](matrices.md#representing-data-in-mlpack). + - If `strict` is `true`, the size of `a` cannot be changed. + - `a` should be a cube type (e.g. `arma::cube`, `arma::fcube`), and `memptr` + should be a non-const pointer of the matrix's element type (e.g. `double*`, + `float*`). + +--- + + * `ClearAlias(a)` + - If `a` is an alias, reset `a` to an empty matrix, without modifying the + aliased memory. `a` is no longer an alias after this call. + +--- + + * `UnwrapAlias(a, in)` + - If `in` is a matrix type (e.g. `arma::mat`), make `a` into an alias of + `in`. + - If `in` is not a matrix type, but instead, e.g., an Armadillo expression, + fill `a` with the results of the evaluated expression `in`. + - This can be used in place of, e.g., `a = in`, to avoid a copy when + possible. + - `a` should be a matrix type that matches the type of the expression or + matrix `in`. + +--- + ### `Range` The `Range` class represents a simple mathematical range (i.e. `[0, 3]`), diff --git a/src/mlpack/core/math/make_alias.hpp b/src/mlpack/core/math/make_alias.hpp index f73578e6ed..93bcd0138a 100644 --- a/src/mlpack/core/math/make_alias.hpp +++ b/src/mlpack/core/math/make_alias.hpp @@ -16,89 +16,69 @@ namespace mlpack { /** - * Make an alias of a dense cube. If strict is true, then the alias cannot be - * resized or pointed at new memory. + * Reconstruct `m` as an alias around the memory `newMem`, with size `numRows` x + * `numCols`. */ -template -arma::Cube MakeAlias(arma::Cube& input, - const bool strict = true) +template +void MakeAlias(MatType& m, + typename MatType::elem_type* newMem, + const size_t numRows, + const size_t numCols, + const bool strict = true, + const typename std::enable_if_t::value>* = 0) { - // Use the advanced constructor. - return arma::Cube(input.memptr(), input.n_rows, input.n_cols, - input.n_slices, false, strict); + // We use placement new to reinitialize the object, since the copy and move + // assignment operators in Armadillo will end up copying memory instead of + // making an alias. + m.~MatType(); + new (&m) MatType(newMem, numRows, numCols, false, strict); } /** - * Make an alias of a dense matrix. If strict is true, then the alias cannot be - * resized or pointed at new memory. + * Reconstruct `c` as an alias around the memory` newMem`, with size `numRows` x + * `numCols` x `numSlices`. */ -template -arma::Mat MakeAlias(arma::Mat& input, - const bool strict = true) +template +void MakeAlias(CubeType& c, + typename CubeType::elem_type* newMem, + const size_t numRows, + const size_t numCols, + const size_t numSlices, + const bool strict = true, + const typename std::enable_if_t::value>* = 0) { - // Use the advanced constructor. - return arma::Mat(input.memptr(), input.n_rows, input.n_cols, false, - strict); + // We use placement new to reinitialize the object, since the copy and move + // assignment operators in Armadillo will end up copying memory instead of + // making an alias. + c.~CubeType(); + new (&c) CubeType(newMem, numRows, numCols, numSlices, false, strict); } /** - * Make an alias of a dense row. If strict is true, then the alias cannot be - * resized or pointed at new memory. + * Make `m` an alias of `in`, using the given size. */ -template -arma::Row MakeAlias(arma::Row& input, - const bool strict = true) +template +void MakeAlias(arma::Mat& m, + const arma::Mat& in, + const size_t numRows, + const size_t numCols, + const bool strict = true) { - // Use the advanced constructor. - return arma::Row(input.memptr(), input.n_elem, false, strict); + MakeAlias(m, (eT*) in.memptr(), numRows, numCols, strict); } /** - * Make an alias of a dense column. If strict is true, then the alias cannot be - * resized or pointed at new memory. + * Make `m` an alias of `in`, using the given size. */ -template -arma::Col MakeAlias(arma::Col& input, - const bool strict = true) +template +void MakeAlias(arma::SpMat& m, + const arma::SpMat& in, + const size_t /* numRows */, + const size_t /* numCols */, + const bool /* strict */) { - // Use the advanced constructor. - return arma::Col(input.memptr(), input.n_elem, false, strict); -} - -/** - * Make a copy of a sparse matrix (an alias is not possible). The strict - * parameter is ignored. - */ -template -arma::SpMat MakeAlias(const arma::SpMat& input, - const bool /* strict */ = true) -{ - // Make a copy... - return arma::SpMat(input); -} - -/** - * Make a copy of a sparse row (an alias is not possible). The strict - * parameter is ignored. - */ -template -arma::SpRow MakeAlias(const arma::SpRow& input, - const bool /* strict */ = true) -{ - // Make a copy... - return arma::SpRow(input); -} - -/** - * Make a copy of a sparse column (an alias is not possible). The strict - * parameter is ignored. - */ -template -arma::SpCol MakeAlias(const arma::SpCol& input, - const bool /* strict */ = true) -{ - // Make a copy... - return arma::SpCol(input); + // We can't make aliases of sparse objects, so just copy it. + m = in; } /** @@ -113,16 +93,15 @@ void ClearAlias(arma::Mat& mat) } /** - * Clear an alias for a sparse matrix. This does nothing because no sparse - * matrices can have aliases. + * Clear an alias so that no data is overwritten. This resets the matrix if it + * is an alias (and does nothing otherwise). */ template void ClearAlias(arma::SpMat& /* mat */) { - // Nothing to do. + // We cannot make aliases of sparse matrices, so, nothing to do. } - } // namespace mlpack #endif diff --git a/src/mlpack/core/math/math.hpp b/src/mlpack/core/math/math.hpp index 66fcfbdd59..3e254863fa 100644 --- a/src/mlpack/core/math/math.hpp +++ b/src/mlpack/core/math/math.hpp @@ -25,5 +25,6 @@ #include "range.hpp" #include "shuffle_data.hpp" #include "trigamma.hpp" +#include "unwrap_alias.hpp" #endif diff --git a/src/mlpack/core/math/unwrap_alias.hpp b/src/mlpack/core/math/unwrap_alias.hpp new file mode 100644 index 0000000000..7afaf30d96 --- /dev/null +++ b/src/mlpack/core/math/unwrap_alias.hpp @@ -0,0 +1,38 @@ +/** + * @file core/math/unwrap_alias.hpp + * @author Ryan Curtin + * + * Make an alias of a matrix if possible, or unwrap it if an expression is + * given. + * + * mlpack is free software; you may redistribute it and/or modify it under the + * terms of the 3-clause BSD license. You should have received a copy of the + * 3-clause BSD license along with mlpack. If not, see + * http://www.opensource.org/licenses/BSD-3-Clause for more information. + */ +#ifndef MLPACK_CORE_MATH_UNWRAP_ALIAS_HPP +#define MLPACK_CORE_MATH_UNWRAP_ALIAS_HPP + +namespace mlpack { + +/** + * If `in` is an expression, unwrap it into `m`. If `in` is a matrix, then + * create `m` as an alias of `in`. + * `numCols`. + */ +template +void UnwrapAlias(MatType& m, const MatType& in) +{ + MakeAlias(m, in, in.n_rows, in.n_cols); +} + +template +void UnwrapAlias(MatType& m, + const InMatType& in) +{ + m = in; +} + +} // namespace mlpack + +#endif diff --git a/src/mlpack/methods/ann/ann.hpp b/src/mlpack/methods/ann/ann.hpp index 9c20ec46a3..e2e2a468d3 100644 --- a/src/mlpack/methods/ann/ann.hpp +++ b/src/mlpack/methods/ann/ann.hpp @@ -17,7 +17,6 @@ #define MLPACK_METHODS_ANN_ANN_HPP #include "forward_decls.hpp" -#include "make_alias.hpp" #include "activation_functions/activation_functions.hpp" #include "augmented/augmented.hpp" diff --git a/src/mlpack/methods/ann/ffn_impl.hpp b/src/mlpack/methods/ann/ffn_impl.hpp index aaf7545606..f95eec3ea1 100644 --- a/src/mlpack/methods/ann/ffn_impl.hpp +++ b/src/mlpack/methods/ann/ffn_impl.hpp @@ -15,8 +15,6 @@ // In case it hasn't been included yet. #include "ffn.hpp" -#include "make_alias.hpp" - namespace mlpack { template& InputDimensions() const { return inputDimensions; } diff --git a/src/mlpack/methods/ann/layer/multi_layer.hpp b/src/mlpack/methods/ann/layer/multi_layer.hpp index 991bafa31b..bad1ef6038 100644 --- a/src/mlpack/methods/ann/layer/multi_layer.hpp +++ b/src/mlpack/methods/ann/layer/multi_layer.hpp @@ -12,7 +12,6 @@ #ifndef MLPACK_METHODS_ANN_LAYER_MULTI_LAYER_HPP #define MLPACK_METHODS_ANN_LAYER_MULTI_LAYER_HPP -#include "../make_alias.hpp" #include "layer.hpp" namespace mlpack { diff --git a/src/mlpack/methods/ann/make_alias.hpp b/src/mlpack/methods/ann/make_alias.hpp deleted file mode 100644 index a0de66a004..0000000000 --- a/src/mlpack/methods/ann/make_alias.hpp +++ /dev/null @@ -1,58 +0,0 @@ -/** - * @file make_alias.hpp - * @author Ryan Curtin - * - * Implementation of `MakeAlias()`, a utility function. This is meant to be - * used in `SetWeights()` calls in various layers, to wrap internal weight - * objects as aliases around the given memory pointers. - * - * mlpack is free software; you may redistribute it and/or modify it under the - * terms of the 3-clause BSD license. You should have received a copy of the - * 3-clause BSD license along with mlpack. If not, see - * http://www.opensource.org/licenses/BSD-3-Clause for more information. - */ -#ifndef MLPACK_METHODS_ANN_MAKE_ALIAS_HPP -#define MLPACK_METHODS_ANN_MAKE_ALIAS_HPP - -#include - -namespace mlpack { - -/** - * Reconstruct `m` as an alias around the memory `newMem`, with size `numRows` x - * `numCols`. - */ -template -void MakeAlias(MatType& m, - typename MatType::elem_type* newMem, - const size_t numRows, - const size_t numCols) -{ - // We use placement new to reinitialize the object, since the copy and move - // assignment operators in Armadillo will end up copying memory instead of - // making an alias. - m.~MatType(); - new (&m) MatType(newMem, numRows, numCols, false, true); -} - -/** - * Reconstruct `c` as an alias around the memory` newMem`, with size `numRows` x - * `numCols` x `numSlices`. - */ -template -void MakeAlias(CubeType& c, - typename CubeType::elem_type* newMem, - const size_t numRows, - const size_t numCols, - const size_t numSlices) -{ - // We use placement new to reinitialize the object, since the copy and move - // assignment operators in Armadillo will end up copying memory instead of - // making an alias. - c.~CubeType(); - new (&c) CubeType(newMem, numRows, numCols, numSlices, false, true); -} - -} // namespace mlpack - -#endif diff --git a/src/mlpack/methods/bias_svd/bias_svd_function_impl.hpp b/src/mlpack/methods/bias_svd/bias_svd_function_impl.hpp index 729b3f3c98..735b2ec66a 100644 --- a/src/mlpack/methods/bias_svd/bias_svd_function_impl.hpp +++ b/src/mlpack/methods/bias_svd/bias_svd_function_impl.hpp @@ -19,13 +19,14 @@ namespace mlpack { template -BiasSVDFunction::BiasSVDFunction(const MatType& data, +BiasSVDFunction::BiasSVDFunction(const MatType& dataIn, const size_t rank, const double lambda) : - data(MakeAlias(const_cast(data), false)), rank(rank), lambda(lambda) { + MakeAlias(data, dataIn, dataIn.n_rows, dataIn.n_cols, false); + // Number of users and items in the data. numUsers = max(data.row(0)) + 1; numItems = max(data.row(1)) + 1; diff --git a/src/mlpack/methods/linear_svm/linear_svm_function_impl.hpp b/src/mlpack/methods/linear_svm/linear_svm_function_impl.hpp index e403f060a2..ea3e8d56ff 100644 --- a/src/mlpack/methods/linear_svm/linear_svm_function_impl.hpp +++ b/src/mlpack/methods/linear_svm/linear_svm_function_impl.hpp @@ -24,18 +24,19 @@ namespace mlpack { template LinearSVMFunction::LinearSVMFunction( - const MatType& dataset, + const MatType& datasetIn, const arma::Row& labels, const size_t numClasses, const double lambda, const double delta, const bool fitIntercept) : - dataset(MakeAlias(const_cast(dataset), false)), numClasses(numClasses), lambda(lambda), delta(delta), fitIntercept(fitIntercept) { + MakeAlias(dataset, datasetIn, datasetIn.n_rows, datasetIn.n_cols, false); + InitializeWeights(initialPoint, dataset.n_rows, numClasses, fitIntercept); initialPoint *= 0.005; diff --git a/src/mlpack/methods/lmnn/lmnn_function_impl.hpp b/src/mlpack/methods/lmnn/lmnn_function_impl.hpp index c2ec0dc7ff..b1779f1731 100644 --- a/src/mlpack/methods/lmnn/lmnn_function_impl.hpp +++ b/src/mlpack/methods/lmnn/lmnn_function_impl.hpp @@ -19,23 +19,24 @@ namespace mlpack { template -LMNNFunction::LMNNFunction(const arma::mat& dataset, - const arma::Row& labels, +LMNNFunction::LMNNFunction(const arma::mat& datasetIn, + const arma::Row& labelsIn, size_t k, double regularization, size_t range, MetricType metric) : - dataset(MakeAlias(const_cast(dataset), false)), - labels(MakeAlias(const_cast&>(labels), false)), k(k), metric(metric), regularization(regularization), iteration(0), range(range), - constraint(dataset, labels, k), - points(dataset.n_cols), + constraint(datasetIn, labelsIn, k), + points(datasetIn.n_cols), impBounds(false) { + MakeAlias(dataset, datasetIn, datasetIn.n_rows, datasetIn.n_cols, false); + MakeAlias(labels, labelsIn, labelsIn.n_rows, labelsIn.n_cols, false); + // Initialize the initial learning point. initialPoint.eye(dataset.n_rows, dataset.n_rows); // Initialize transformed dataset to base dataset. diff --git a/src/mlpack/methods/logistic_regression/logistic_regression_function_impl.hpp b/src/mlpack/methods/logistic_regression/logistic_regression_function_impl.hpp index 4e861ac287..3b7c2282b5 100644 --- a/src/mlpack/methods/logistic_regression/logistic_regression_function_impl.hpp +++ b/src/mlpack/methods/logistic_regression/logistic_regression_function_impl.hpp @@ -21,15 +21,17 @@ namespace mlpack { template LogisticRegressionFunction::LogisticRegressionFunction( - const MatType& predictors, - const arma::Row& responses, + const MatType& predictorsIn, + const arma::Row& responsesIn, const double lambda) : - // We promise to be well-behaved... the elements won't be modified. - predictors(MakeAlias(const_cast(predictors), false)), - responses(MakeAlias(const_cast&>(responses), - false)), lambda(lambda) { + // We promise to be well-behaved... the elements won't be modified. + MakeAlias(this->predictors, predictorsIn, predictorsIn.n_rows, + predictorsIn.n_cols, false); + MakeAlias(this->responses, responsesIn, responsesIn.n_rows, + responsesIn.n_cols, false); + // Sanity check. if (responses.n_elem != predictors.n_cols) { diff --git a/src/mlpack/methods/nca/nca_softmax_error_function_impl.hpp b/src/mlpack/methods/nca/nca_softmax_error_function_impl.hpp index c83f3afb7f..1ce40fb853 100644 --- a/src/mlpack/methods/nca/nca_softmax_error_function_impl.hpp +++ b/src/mlpack/methods/nca/nca_softmax_error_function_impl.hpp @@ -22,14 +22,15 @@ namespace mlpack { // Initialize with the given kernel. template SoftmaxErrorFunction::SoftmaxErrorFunction( - const arma::mat& dataset, - const arma::Row& labels, + const arma::mat& datasetIn, + const arma::Row& labelsIn, MetricType metric) : - dataset(MakeAlias(const_cast(dataset), false)), - labels(MakeAlias(const_cast&>(labels), false)), metric(metric), precalculated(false) -{ /* nothing to do */ } +{ + MakeAlias(dataset, datasetIn, datasetIn.n_rows, datasetIn.n_cols, false); + MakeAlias(labels, labelsIn, labelsIn.n_rows, labelsIn.n_cols, false); +} //! Shuffle the dataset. template diff --git a/src/mlpack/methods/regularized_svd/regularized_svd_function_impl.hpp b/src/mlpack/methods/regularized_svd/regularized_svd_function_impl.hpp index a304e14598..cdff2ae422 100644 --- a/src/mlpack/methods/regularized_svd/regularized_svd_function_impl.hpp +++ b/src/mlpack/methods/regularized_svd/regularized_svd_function_impl.hpp @@ -18,13 +18,14 @@ namespace mlpack { template -RegularizedSVDFunction::RegularizedSVDFunction(const MatType& data, +RegularizedSVDFunction::RegularizedSVDFunction(const MatType& dataIn, const size_t rank, const double lambda) : - data(MakeAlias(const_cast(data), false)), rank(rank), lambda(lambda) { + MakeAlias(data, dataIn, dataIn.n_rows, dataIn.n_cols, false); + // Number of users and items in the data. numUsers = max(data.row(0)) + 1; numItems = max(data.row(1)) + 1; diff --git a/src/mlpack/methods/softmax_regression/softmax_regression_function_impl.hpp b/src/mlpack/methods/softmax_regression/softmax_regression_function_impl.hpp index 70137a079f..2e0da0ecda 100644 --- a/src/mlpack/methods/softmax_regression/softmax_regression_function_impl.hpp +++ b/src/mlpack/methods/softmax_regression/softmax_regression_function_impl.hpp @@ -18,16 +18,17 @@ namespace mlpack { template inline SoftmaxRegressionFunction::SoftmaxRegressionFunction( - const MatType& data, + const MatType& dataIn, const arma::Row& labels, const size_t numClasses, const double lambda, const bool fitIntercept) : - data(MakeAlias(const_cast(data), false)), numClasses(numClasses), lambda(lambda), fitIntercept(fitIntercept) { + MakeAlias(data, dataIn, dataIn.n_rows, dataIn.n_cols, false); + // Initialize the parameters to suitable values. initialPoint = InitializeWeights(); diff --git a/src/mlpack/methods/svdplusplus/svdplusplus_function_impl.hpp b/src/mlpack/methods/svdplusplus/svdplusplus_function_impl.hpp index 08b7ac9952..2b6eaa3e0e 100644 --- a/src/mlpack/methods/svdplusplus/svdplusplus_function_impl.hpp +++ b/src/mlpack/methods/svdplusplus/svdplusplus_function_impl.hpp @@ -20,15 +20,16 @@ namespace mlpack { template SVDPlusPlusFunction::SVDPlusPlusFunction( - const MatType& data, + const MatType& dataIn, const arma::sp_mat& implicitData, const size_t rank, const double lambda) : - data(MakeAlias(const_cast(data), false)), implicitData(implicitData), rank(rank), lambda(lambda) { + MakeAlias(data, dataIn, dataIn.n_rows, dataIn.n_cols, false); + // Number of users and items in the data. numUsers = max(data.row(0)) + 1; numItems = max(data.row(1)) + 1; From 41ad950d5e13b4a794782aa82973b8f4465a8436 Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Tue, 2 Apr 2024 14:52:15 -0400 Subject: [PATCH 02/11] Templatize PCA and its dependencies to allow different matrix types. --- .../core/tree/cosine_tree/cosine_tree.hpp | 61 ++++---- .../tree/cosine_tree/cosine_tree_impl.hpp | 125 +++++++++------- src/mlpack/core/util/arma_traits.hpp | 45 ++++++ .../randomized_block_krylov_svd.hpp | 18 ++- .../randomized_block_krylov_svd_impl.hpp | 34 +++-- .../quic_svd_method.hpp | 2 +- .../exact_svd_method.hpp | 17 ++- .../quic_svd_method.hpp | 17 ++- .../randomized_block_krylov_method.hpp | 15 +- .../randomized_svd_method.hpp | 15 +- src/mlpack/methods/pca/pca.hpp | 39 +++-- src/mlpack/methods/pca/pca_impl.hpp | 132 ++++++++++++----- src/mlpack/methods/quic_svd/quic_svd.hpp | 27 ++-- src/mlpack/methods/quic_svd/quic_svd_impl.hpp | 50 ++++--- .../methods/randomized_svd/randomized_svd.hpp | 42 +++--- .../randomized_svd/randomized_svd_impl.hpp | 61 ++++---- src/mlpack/tests/cosine_tree_test.cpp | 50 +++---- src/mlpack/tests/pca_test.cpp | 137 ++++++++++++++++++ src/mlpack/tests/policy_gradient_test.cpp | 4 +- src/mlpack/tests/quic_svd_test.cpp | 6 +- 20 files changed, 610 insertions(+), 287 deletions(-) diff --git a/src/mlpack/core/tree/cosine_tree/cosine_tree.hpp b/src/mlpack/core/tree/cosine_tree/cosine_tree.hpp index 5966688e37..cb5b52fe25 100644 --- a/src/mlpack/core/tree/cosine_tree/cosine_tree.hpp +++ b/src/mlpack/core/tree/cosine_tree/cosine_tree.hpp @@ -20,14 +20,20 @@ namespace mlpack { // Predeclare classes for CosineNodeQueue typedef. class CompareCosineNode; + +template class CosineTree; // CosineNodeQueue typedef. -typedef std::vector CosineNodeQueue; +template +using CosineNodeQueue = std::vector*>; +template class CosineTree { public: + typedef typename GetDenseColType::type VecType; + /** * CosineTree constructor for the root node of the tree. It initializes the * necessary variables required for splitting of the node, and building the @@ -36,7 +42,7 @@ class CosineTree * * @param dataset Matrix for which cosine tree is constructed. */ - CosineTree(const arma::mat& dataset); + CosineTree(const MatType& dataset); /** * CosineTree constructor for nodes other than the root node of the tree. It @@ -63,7 +69,7 @@ class CosineTree * @param epsilon Error tolerance fraction for calculated subspace. * @param delta Cumulative probability for Monte Carlo error lower bound. */ - CosineTree(const arma::mat& dataset, + CosineTree(const MatType& dataset, const double epsilon, const double delta); @@ -110,10 +116,10 @@ class CosineTree * @param newBasisVector Orthonormalized centroid of the node. * @param addBasisVector Address to additional basis vector. */ - void ModifiedGramSchmidt(CosineNodeQueue& treeQueue, - arma::vec& centroid, - arma::vec& newBasisVector, - arma::vec* addBasisVector = NULL); + void ModifiedGramSchmidt(CosineNodeQueue& treeQueue, + VecType& centroid, + VecType& newBasisVector, + VecType* addBasisVector = NULL); /** * Estimates the squared error of the projection of the input node's matrix @@ -128,16 +134,16 @@ class CosineTree * @param addBasisVector2 Address to second additional basis vector. */ double MonteCarloError(CosineTree* node, - CosineNodeQueue& treeQueue, - arma::vec* addBasisVector1 = NULL, - arma::vec* addBasisVector2 = NULL); + CosineNodeQueue& treeQueue, + VecType* addBasisVector1 = NULL, + VecType* addBasisVector2 = NULL); /** * Constructs the final basis matrix, after the cosine tree construction. * * @param treeQueue Priority queue of cosine nodes. */ - void ConstructBasis(CosineNodeQueue& treeQueue); + void ConstructBasis(CosineNodeQueue& treeQueue); /** * This function splits the cosine node into two children based on the cosines @@ -153,7 +159,8 @@ class CosineTree * randomly generated values in the range [0, 1]. */ void ColumnSamplesLS(std::vector& sampledIndices, - arma::vec& probabilities, size_t numSamples); + VecType& probabilities, + size_t numSamples); /** * Sample a point from the Length-Squared distribution of the cosine node. The @@ -175,7 +182,9 @@ class CosineTree * @param start Starting index of the distribution interval to search in. * @param end Ending index of the distribution interval to search in. */ - size_t BinarySearch(arma::vec& cDistribution, double value, size_t start, + size_t BinarySearch(VecType& cDistribution, + double value, + size_t start, size_t end); /** @@ -185,7 +194,7 @@ class CosineTree * * @param cosines Vector to store the cosine values in. */ - void CalculateCosines(arma::vec& cosines); + void CalculateCosines(VecType& cosines); /** * Calculate centroid of the columns present in the node. The calculated @@ -194,10 +203,10 @@ class CosineTree void CalculateCentroid(); //! Returns the basis of the constructed subspace. - void GetFinalBasis(arma::mat& finalBasis) { finalBasis = basis; } + void GetFinalBasis(MatType& finalBasis) { finalBasis = basis; } //! Get pointer to the dataset matrix. - const arma::mat& GetDataset() const { return *dataset; } + const MatType& GetDataset() const { return *dataset; } //! Get the indices of columns in the node. std::vector& VectorIndices() { return indices; } @@ -208,13 +217,13 @@ class CosineTree double L2Error() const { return l2Error; } //! Get pointer to the centroid vector. - arma::vec& Centroid() { return centroid; } + VecType& Centroid() { return centroid; } //! Set the basis vector of the node. - void BasisVector(arma::vec& bVector) { this->basisVector = bVector; } + void BasisVector(VecType& bVector) { this->basisVector = bVector; } //! Get the basis vector of the node. - arma::vec& BasisVector() { return basisVector; } + VecType& BasisVector() { return basisVector; } //! Get pointer to the parent node. CosineTree* Parent() const { return parent; } @@ -242,11 +251,11 @@ class CosineTree private: //! Matrix for which cosine tree is constructed. - const arma::mat* dataset; + const MatType* dataset; //! Cumulative probability for Monte Carlo error lower bound. double delta; //! Subspace basis of the input dataset. - arma::mat basis; + MatType basis; //! Parent of the node. CosineTree* parent; //! Left child of the node. @@ -256,11 +265,11 @@ class CosineTree //! Indices of columns of input matrix in the node. std::vector indices; //! L2-norm squared of columns in the node. - arma::vec l2NormsSquared; + VecType l2NormsSquared; //! Centroid of columns of input matrix in the node. - arma::vec centroid; + VecType centroid; //! Orthonormalized basis vector of the node. - arma::vec basisVector; + VecType basisVector; //! Index of split point of cosine node. size_t splitPointIndex; //! Number of columns of input matrix in the node. @@ -277,7 +286,9 @@ class CompareCosineNode { public: // Comparison function for construction of priority queue. - bool operator() (const CosineTree* a, const CosineTree* b) const + template + bool operator() (const CosineTree* a, + const CosineTree* b) const { return a->L2Error() < b->L2Error(); } diff --git a/src/mlpack/core/tree/cosine_tree/cosine_tree_impl.hpp b/src/mlpack/core/tree/cosine_tree/cosine_tree_impl.hpp index dcc036465f..7ef0393399 100644 --- a/src/mlpack/core/tree/cosine_tree/cosine_tree_impl.hpp +++ b/src/mlpack/core/tree/cosine_tree/cosine_tree_impl.hpp @@ -16,7 +16,8 @@ namespace mlpack { -inline CosineTree::CosineTree(const arma::mat& dataset) : +template +inline CosineTree::CosineTree(const MatType& dataset) : dataset(&dataset), parent(NULL), left(NULL), @@ -32,7 +33,7 @@ inline CosineTree::CosineTree(const arma::mat& dataset) : for (size_t i = 0; i < numColumns; ++i) { indices[i] = i; - double l2Norm = arma::norm(dataset.col(i), 2); + double l2Norm = (double) arma::norm(dataset.col(i), 2); l2NormsSquared(i) = l2Norm * l2Norm; } @@ -45,8 +46,9 @@ inline CosineTree::CosineTree(const arma::mat& dataset) : splitPointIndex = ColumnSampleLS(); } -inline CosineTree::CosineTree(CosineTree& parentNode, - const std::vector& subIndices) : +template +inline CosineTree::CosineTree(CosineTree& parentNode, + const std::vector& subIndices) : dataset(&parentNode.GetDataset()), parent(&parentNode), left(NULL), @@ -74,9 +76,10 @@ inline CosineTree::CosineTree(CosineTree& parentNode, splitPointIndex = ColumnSampleLS(); } -inline CosineTree::CosineTree(const arma::mat& dataset, - const double epsilon, - const double delta) : +template +inline CosineTree::CosineTree(const MatType& dataset, + const double epsilon, + const double delta) : dataset(&dataset), delta(delta), left(NULL), @@ -84,15 +87,16 @@ inline CosineTree::CosineTree(const arma::mat& dataset, localDataset(false) { // Declare the cosine tree priority queue. - CosineNodeQueue treeQueue; + CosineNodeQueue treeQueue; CompareCosineNode comp; // Define root node of the tree and add it to the queue. CosineTree root(dataset); - arma::vec tempVector = arma::zeros(dataset.n_rows); + VecType tempVector = arma::zeros(dataset.n_rows); root.L2Error(-1.0); // We don't know what the error is. root.BasisVector(tempVector); - treeQueue.push_back(&root); // treeQueue is empty now, so we don't need to call std::push_heap here. + treeQueue.push_back(&root); + // treeQueue is empty now, so we don't need to call std::push_heap here. // Initialize Monte Carlo error estimate for comparison. double monteCarloError = root.FrobNormSquared(); @@ -128,7 +132,7 @@ inline CosineTree::CosineTree(const arma::mat& dataset, currentRight = currentNode->Right(); // Calculate basis vectors of left and right children. - arma::vec lBasisVector, rBasisVector; + VecType lBasisVector, rBasisVector; ModifiedGramSchmidt(treeQueue, currentLeft->Centroid(), lBasisVector); ModifiedGramSchmidt(treeQueue, currentRight->Centroid(), rBasisVector, @@ -157,9 +161,10 @@ inline CosineTree::CosineTree(const arma::mat& dataset, } //! Copy the given tree. -inline CosineTree::CosineTree(const CosineTree& other) : +template +inline CosineTree::CosineTree(const CosineTree& other) : // Copy matrix, but only if we are the root. - dataset((other.parent == NULL) ? new arma::mat(*other.dataset) : NULL), + dataset((other.parent == NULL) ? new MatType(*other.dataset) : NULL), delta(other.delta), parent(NULL), left(NULL), @@ -210,7 +215,9 @@ inline CosineTree::CosineTree(const CosineTree& other) : } //! Copy assignment operator: copy the given other tree. -inline CosineTree& CosineTree::operator=(const CosineTree& other) +template +inline CosineTree& CosineTree::operator=( + const CosineTree& other) { // Return if it's the same tree. if (this == &other) @@ -224,7 +231,7 @@ inline CosineTree& CosineTree::operator=(const CosineTree& other) delete right; // Performing a deep copy of the dataset. - dataset = (other.parent == NULL) ? new arma::mat(*other.dataset) : NULL; + dataset = (other.parent == NULL) ? new MatType(*other.dataset) : NULL; delta = other.delta; parent = other.Parent(); @@ -278,7 +285,8 @@ inline CosineTree& CosineTree::operator=(const CosineTree& other) } //! Move the given tree. -inline CosineTree::CosineTree(CosineTree&& other) : +template +inline CosineTree::CosineTree(CosineTree&& other) : dataset(other.dataset), delta(std::move(other.delta)), parent(other.parent), @@ -313,7 +321,8 @@ inline CosineTree::CosineTree(CosineTree&& other) : } //! Move assignment operator: take ownership of the given tree. -inline CosineTree& CosineTree::operator=(CosineTree&& other) +template +inline CosineTree& CosineTree::operator=(CosineTree&& other) { // Return if it's the same tree. if (this == &other) @@ -360,7 +369,8 @@ inline CosineTree& CosineTree::operator=(CosineTree&& other) return *this; } -inline CosineTree::~CosineTree() +template +inline CosineTree::~CosineTree() { if (localDataset) delete dataset; @@ -370,17 +380,19 @@ inline CosineTree::~CosineTree() delete right; } -inline void CosineTree::ModifiedGramSchmidt(CosineNodeQueue& treeQueue, - arma::vec& centroid, - arma::vec& newBasisVector, - arma::vec* addBasisVector) +template +inline void CosineTree::ModifiedGramSchmidt( + CosineNodeQueue& treeQueue, + typename CosineTree::VecType& centroid, + typename CosineTree::VecType& newBasisVector, + typename CosineTree::VecType* addBasisVector) { // Set new basis vector to centroid. newBasisVector = centroid; // Variables for iterating throught the priority queue. CosineTree *currentNode; - CosineNodeQueue::const_iterator i = treeQueue.cbegin(); + typename CosineNodeQueue::const_iterator i = treeQueue.cbegin(); // For every vector in the current basis, remove its projection from the // centroid. @@ -388,14 +400,14 @@ inline void CosineTree::ModifiedGramSchmidt(CosineNodeQueue& treeQueue, { currentNode = *i; - double projection = dot(currentNode->BasisVector(), centroid); + double projection = (double) dot(currentNode->BasisVector(), centroid); newBasisVector -= projection * currentNode->BasisVector(); } // If additional basis vector is passed, take it into account. if (addBasisVector) { - double projection = dot(*addBasisVector, centroid); + double projection = (double) dot(*addBasisVector, centroid); newBasisVector -= *addBasisVector * projection; } @@ -404,13 +416,15 @@ inline void CosineTree::ModifiedGramSchmidt(CosineNodeQueue& treeQueue, newBasisVector /= arma::norm(newBasisVector, 2); } -inline double CosineTree::MonteCarloError(CosineTree* node, - CosineNodeQueue& treeQueue, - arma::vec* addBasisVector1, - arma::vec* addBasisVector2) +template +inline double CosineTree::MonteCarloError( + CosineTree* node, + CosineNodeQueue& treeQueue, + typename CosineTree::VecType* addBasisVector1, + typename CosineTree::VecType* addBasisVector2) { std::vector sampledIndices; - arma::vec probabilities; + VecType probabilities; // Sample O(log m) points from the input node's distribution. // 'm' is the number of columns present in the node. @@ -418,10 +432,10 @@ inline double CosineTree::MonteCarloError(CosineTree* node, node->ColumnSamplesLS(sampledIndices, probabilities, numSamples); // Get pointer to the original dataset. - const arma::mat& dataset = node->GetDataset(); + const MatType& dataset = node->GetDataset(); // Initialize weighted projection magnitudes as zeros. - arma::vec weightedMagnitudes; + VecType weightedMagnitudes; weightedMagnitudes.zeros(numSamples); // Set size of projection vector, depending on whether additional basis @@ -436,11 +450,11 @@ inline double CosineTree::MonteCarloError(CosineTree* node, for (size_t i = 0; i < numSamples; ++i) { // Initialize projection as a vector of zeros. - arma::vec projection; + VecType projection; projection.zeros(projectionSize); CosineTree *currentNode; - CosineNodeQueue::const_iterator j = treeQueue.cbegin(); + typename CosineNodeQueue::const_iterator j = treeQueue.cbegin(); size_t k = 0; // Compute the projection of the sampled vector onto the existing subspace. @@ -488,14 +502,16 @@ inline double CosineTree::MonteCarloError(CosineTree* node, return (node->FrobNormSquared() - lowerBound); } -inline void CosineTree::ConstructBasis(CosineNodeQueue& treeQueue) +template +inline void CosineTree::ConstructBasis( + CosineNodeQueue& treeQueue) { // Initialize basis as matrix of zeros. basis.zeros(dataset->n_rows, treeQueue.size()); // Variables for iterating through the priority queue. CosineTree *currentNode; - CosineNodeQueue::const_iterator i = treeQueue.cbegin(); + typename CosineNodeQueue::const_iterator i = treeQueue.cbegin(); // Transfer basis vectors from the queue to the basis matrix. size_t j = 0; @@ -506,7 +522,8 @@ inline void CosineTree::ConstructBasis(CosineNodeQueue& treeQueue) } } -inline void CosineTree::CosineNodeSplit() +template +inline void CosineTree::CosineNodeSplit() { // If less than two points, splitting does not make sense---there is nothing // to split. @@ -514,7 +531,7 @@ inline void CosineTree::CosineNodeSplit() return; // Calculate cosines with respect to the splitting point. - arma::vec cosines; + VecType cosines; CalculateCosines(cosines); // Compute maximum and minimum cosine values. @@ -543,12 +560,14 @@ inline void CosineTree::CosineNodeSplit() right = new CosineTree(*this, rightIndices); } -inline void CosineTree::ColumnSamplesLS(std::vector& sampledIndices, - arma::vec& probabilities, - size_t numSamples) +template +inline void CosineTree::ColumnSamplesLS( + std::vector& sampledIndices, + typename CosineTree::VecType& probabilities, + size_t numSamples) { // Initialize the cumulative distribution vector size. - arma::vec cDistribution; + VecType cDistribution; cDistribution.zeros(numColumns + 1); // Calculate cumulative length-squared distribution for the node. @@ -575,7 +594,8 @@ inline void CosineTree::ColumnSamplesLS(std::vector& sampledIndices, } } -inline size_t CosineTree::ColumnSampleLS() +template +inline size_t CosineTree::ColumnSampleLS() { // If only one element is present, there can only be one sample. if (numColumns < 2) @@ -584,7 +604,7 @@ inline size_t CosineTree::ColumnSampleLS() } // Initialize the cumulative distribution vector size. - arma::vec cDistribution; + VecType cDistribution; cDistribution.zeros(numColumns + 1); // Calculate cumulative length-squared distribution for the node. @@ -602,10 +622,12 @@ inline size_t CosineTree::ColumnSampleLS() return BinarySearch(cDistribution, randValue, start, end); } -inline size_t CosineTree::BinarySearch(arma::vec& cDistribution, - double value, - size_t start, - size_t end) +template +inline size_t CosineTree::BinarySearch( + typename CosineTree::VecType& cDistribution, + double value, + size_t start, + size_t end) { size_t pivot = (start + end) / 2; @@ -630,7 +652,9 @@ inline size_t CosineTree::BinarySearch(arma::vec& cDistribution, } } -inline void CosineTree::CalculateCosines(arma::vec& cosines) +template +inline void CosineTree::CalculateCosines( + typename CosineTree::VecType& cosines) { // Initialize cosine vector as a vector of zeros. cosines.zeros(numColumns); @@ -652,7 +676,8 @@ inline void CosineTree::CalculateCosines(arma::vec& cosines) } } -inline void CosineTree::CalculateCentroid() +template +inline void CosineTree::CalculateCentroid() { // Initialize centroid as vector of zeros. centroid.zeros(dataset->n_rows); diff --git a/src/mlpack/core/util/arma_traits.hpp b/src/mlpack/core/util/arma_traits.hpp index 923ea3da72..2dd360761b 100644 --- a/src/mlpack/core/util/arma_traits.hpp +++ b/src/mlpack/core/util/arma_traits.hpp @@ -254,4 +254,49 @@ struct GetSparseMatType> typedef arma::SpMat type; }; +// Get whether or not the given type is a base matrix type (e.g. not an +// expression). + +template +struct IsBaseMatType +{ + constexpr static bool value = false; +}; + +template +struct IsBaseMatType> +{ + constexpr static bool value = true; +}; + +template +struct IsBaseMatType> +{ + constexpr static bool value = true; +}; + +template +struct IsBaseMatType> +{ + constexpr static bool value = true; +}; + +template +struct IsBaseMatType> +{ + constexpr static bool value = true; +}; + +template +struct IsBaseMatType> +{ + constexpr static bool value = true; +}; + +template +struct IsBaseMatType> +{ + constexpr static bool value = true; +}; + #endif diff --git a/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd.hpp b/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd.hpp index add30fb990..ad32145c06 100644 --- a/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd.hpp +++ b/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd.hpp @@ -69,10 +69,11 @@ class RandomizedBlockKrylovSVD * @param rank Rank of the approximation (Default: number of rows.) * @param blockSize The block size, must be >= rank (Default: rank + 10). */ - RandomizedBlockKrylovSVD(const arma::mat& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, + template + RandomizedBlockKrylovSVD(const InMatType& data, + MatType& u, + VecType& s, + MatType& v, const size_t maxIterations = 2, const size_t rank = 0, const size_t blockSize = 0); @@ -97,10 +98,11 @@ class RandomizedBlockKrylovSVD * @param s Diagonal matrix of singular values. * @param rank Rank of the approximation. */ - void Apply(const arma::mat& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, + template + void Apply(const InMatType& data, + MatType& u, + VecType& s, + MatType& v, const size_t rank); //! Get the number of iterations for the power method. diff --git a/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp b/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp index 982321ecbf..a042546252 100644 --- a/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp +++ b/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp @@ -16,11 +16,12 @@ namespace mlpack { +template inline RandomizedBlockKrylovSVD::RandomizedBlockKrylovSVD( - const arma::mat& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, + const InMatType& data, + MatType& u, + VecType& s, + MatType& v, const size_t maxIterations, const size_t rank, const size_t blockSize) : @@ -46,13 +47,14 @@ inline RandomizedBlockKrylovSVD::RandomizedBlockKrylovSVD( /* Nothing to do here */ } -inline void RandomizedBlockKrylovSVD::Apply(const arma::mat& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, +template +inline void RandomizedBlockKrylovSVD::Apply(const InMatType& data, + MatType& u, + VecType& s, + MatType& v, const size_t rank) { - arma::mat Q, R, block, blockIteration; + MatType Q, R, block, blockIteration; if (blockSize == 0) { @@ -60,28 +62,28 @@ inline void RandomizedBlockKrylovSVD::Apply(const arma::mat& data, } // Random block initialization. - arma::mat G = arma::randn(data.n_cols, blockSize); + MatType G = arma::randn(data.n_cols, blockSize); // Construct and orthonormalize Krylov subspace. - arma::mat K(data.n_rows, blockSize * (maxIterations + 1)); + MatType K(data.n_rows, blockSize * (maxIterations + 1)); // Create a working matrix using data from writable auxiliary memory // (K matrix). Doing so avoids an uncessary copy in upcoming step. - block = arma::mat(K.memptr(), data.n_rows, blockSize, false, false); + MakeAlias(block, K.memptr(), data.n_rows, blockSize, false); arma::qr_econ(block, R, data * G); for (size_t blockOffset = block.n_elem; blockOffset < K.n_elem; blockOffset += block.n_elem) { // Temporary working matrix to store the result in the correct place. - blockIteration = arma::mat(K.memptr() + blockOffset, block.n_rows, - block.n_cols, false, false); + MakeAlias(blockIteration, K.memptr() + blockOffset, block.n_rows, + block.n_cols, false); arma::qr_econ(blockIteration, R, data * (data.t() * block)); // Update working matrix for the next iteration. - block = arma::mat(K.memptr() + blockOffset, block.n_rows, block.n_cols, - false, false); + MakeAlias(block, K.memptr() + blockOffset, block.n_rows, block.n_cols, + false); } arma::qr_econ(Q, R, K); diff --git a/src/mlpack/methods/cf/decomposition_policies/quic_svd_method.hpp b/src/mlpack/methods/cf/decomposition_policies/quic_svd_method.hpp index 775a7c8803..a8b6cdbe39 100644 --- a/src/mlpack/methods/cf/decomposition_policies/quic_svd_method.hpp +++ b/src/mlpack/methods/cf/decomposition_policies/quic_svd_method.hpp @@ -74,7 +74,7 @@ class QUIC_SVDPolicy arma::mat data(cleanedData); // Do singular value decomposition using the quic SVD algorithm. - QUIC_SVD quicsvd; + QUIC_SVD<> quicsvd; quicsvd.Apply(data, w, h, sigma); // Sigma matrix is multiplied to w. diff --git a/src/mlpack/methods/pca/decomposition_policies/exact_svd_method.hpp b/src/mlpack/methods/pca/decomposition_policies/exact_svd_method.hpp index 4f42d23eac..cd37fee403 100644 --- a/src/mlpack/methods/pca/decomposition_policies/exact_svd_method.hpp +++ b/src/mlpack/methods/pca/decomposition_policies/exact_svd_method.hpp @@ -37,19 +37,20 @@ class ExactSVDPolicy * @param eigvec Matrix to put eigenvectors (loadings) into. * @param * (rank) Rank of the decomposition. */ - void Apply(const arma::mat& data, - const arma::mat& centeredData, - arma::mat& transformedData, - arma::vec& eigVal, - arma::mat& eigvec, + template + void Apply(const InMatType& /* data */, + const MatType& centeredData, + MatType& transformedData, + VecType& eigVal, + MatType& eigvec, const size_t /* rank */) { // This matrix will store the right singular vectors; we do not need them. - arma::mat v; + MatType v; // Do singular value decomposition. Use the economical singular value // decomposition if the columns are much larger than the rows. - if (data.n_rows < data.n_cols) + if (centeredData.n_rows < centeredData.n_cols) { // Do economical singular value decomposition and compute only the left // singular vectors. @@ -63,7 +64,7 @@ class ExactSVDPolicy // Now we must square the singular values to get the eigenvalues. // In addition we must divide by the number of points, because the // covariance matrix is X * X' / (N - 1). - eigVal %= eigVal / (data.n_cols - 1); + eigVal %= eigVal / (centeredData.n_cols - 1); // Project the samples to the principals. transformedData = trans(eigvec) * centeredData; diff --git a/src/mlpack/methods/pca/decomposition_policies/quic_svd_method.hpp b/src/mlpack/methods/pca/decomposition_policies/quic_svd_method.hpp index 6a84ccbef6..240cbe3e1c 100644 --- a/src/mlpack/methods/pca/decomposition_policies/quic_svd_method.hpp +++ b/src/mlpack/methods/pca/decomposition_policies/quic_svd_method.hpp @@ -49,23 +49,24 @@ class QUICSVDPolicy * @param eigvec Matrix to put eigenvectors (loadings) into. * @param * (rank) Rank of the decomposition. */ - void Apply(const arma::mat& data, - const arma::mat& centeredData, - arma::mat& transformedData, - arma::vec& eigVal, - arma::mat& eigvec, + template + void Apply(const InMatType& /* data */, + const MatType& centeredData, + MatType& transformedData, + VecType& eigVal, + MatType& eigvec, const size_t /* rank */) { // This matrix will store the right singular vectors; we do not need them. - arma::mat v, sigma; + MatType v, sigma; // Do singular value decomposition using the QUIC-SVD algorithm. - QUIC_SVD quicsvd(centeredData, eigvec, v, sigma, epsilon, delta); + QUIC_SVD quicsvd(centeredData, eigvec, v, sigma, epsilon, delta); // Now we must square the singular values to get the eigenvalues. // In addition we must divide by the number of points, because the // covariance matrix is X * X' / (N - 1). - eigVal = pow(arma::diagvec(sigma), 2) / (data.n_cols - 1); + eigVal = pow(arma::diagvec(sigma), 2) / (centeredData.n_cols - 1); // Project the samples to the principals. transformedData = trans(eigvec) * centeredData; diff --git a/src/mlpack/methods/pca/decomposition_policies/randomized_block_krylov_method.hpp b/src/mlpack/methods/pca/decomposition_policies/randomized_block_krylov_method.hpp index 7836b72c8a..64aad7b4b8 100644 --- a/src/mlpack/methods/pca/decomposition_policies/randomized_block_krylov_method.hpp +++ b/src/mlpack/methods/pca/decomposition_policies/randomized_block_krylov_method.hpp @@ -52,15 +52,16 @@ class RandomizedBlockKrylovSVDPolicy * @param eigvec Matrix to put eigenvectors (loadings) into. * @param rank Rank of the decomposition. */ - void Apply(const arma::mat& data, - const arma::mat& centeredData, - arma::mat& transformedData, - arma::vec& eigVal, - arma::mat& eigvec, + template + void Apply(const InMatType& /* data */, + const MatType& centeredData, + MatType& transformedData, + VecType& eigVal, + MatType& eigvec, const size_t rank) { // This matrix will store the right singular vectors; we do not need them. - arma::mat v; + MatType v; // Do singular value decomposition using the randomized block krylov SVD // algorithm. @@ -70,7 +71,7 @@ class RandomizedBlockKrylovSVDPolicy // Now we must square the singular values to get the eigenvalues. // In addition we must divide by the number of points, because the // covariance matrix is X * X' / (N - 1). - eigVal %= eigVal / (data.n_cols - 1); + eigVal %= eigVal / (centeredData.n_cols - 1); // Project the samples to the principals. transformedData = trans(eigvec) * centeredData; diff --git a/src/mlpack/methods/pca/decomposition_policies/randomized_svd_method.hpp b/src/mlpack/methods/pca/decomposition_policies/randomized_svd_method.hpp index ceb2fab0cd..cec751717e 100644 --- a/src/mlpack/methods/pca/decomposition_policies/randomized_svd_method.hpp +++ b/src/mlpack/methods/pca/decomposition_policies/randomized_svd_method.hpp @@ -53,15 +53,16 @@ class RandomizedSVDPCAPolicy * @param eigvec Matrix to put eigenvectors (loadings) into. * @param rank Rank of the decomposition. */ - void Apply(const arma::mat& data, - const arma::mat& centeredData, - arma::mat& transformedData, - arma::vec& eigVal, - arma::mat& eigvec, + template + void Apply(const InMatType& data, + const MatType& centeredData, + MatType& transformedData, + VecType& eigVal, + MatType& eigvec, const size_t rank) { // This matrix will store the right singular vectors; we do not need them. - arma::mat v; + MatType v; // Do singular value decomposition using the randomized SVD algorithm. RandomizedSVD rsvd(iteratedPower, maxIterations); @@ -70,7 +71,7 @@ class RandomizedSVDPCAPolicy // Now we must square the singular values to get the eigenvalues. // In addition we must divide by the number of points, because the // covariance matrix is X * X' / (N - 1). - eigVal %= eigVal / (data.n_cols - 1); + eigVal %= eigVal / (centeredData.n_cols - 1); // Project the samples to the principals. transformedData = trans(eigvec) * centeredData; diff --git a/src/mlpack/methods/pca/pca.hpp b/src/mlpack/methods/pca/pca.hpp index 295b5a188c..612fdebbbf 100644 --- a/src/mlpack/methods/pca/pca.hpp +++ b/src/mlpack/methods/pca/pca.hpp @@ -52,10 +52,13 @@ class PCA * @param eigVal Vector to put eigenvalues into. * @param eigvec Matrix to put eigenvectors (loadings) into. */ - void Apply(const arma::mat& data, - arma::mat& transformedData, - arma::vec& eigVal, - arma::mat& eigvec); + template + void Apply(const MatType& data, + OutMatType& transformedData, + VecType& eigVal, + OutMatType& eigvec); /** * Apply Principal Component Analysis to the provided data set. It is safe @@ -65,17 +68,21 @@ class PCA * @param transformedData Matrix to store results of PCA in. * @param eigVal Vector to put eigenvalues into. */ - void Apply(const arma::mat& data, - arma::mat& transformedData, - arma::vec& eigVal); + template + void Apply(const MatType& data, + OutMatType& transformedData, + VecType& eigVal); /** * Apply Principal Component Analysis to the provided data set. It is safe * to pass the same matrix reference for both data and transformedData. * @param data Data matrix. * @param transformedData Matrix to store results of PCA in. */ - void Apply(const arma::mat& data, - arma::mat& transformedData); + template + void Apply(const MatType& data, + OutMatType& transformedData); /** * Use PCA for dimensionality reduction on the given dataset. This will save @@ -88,10 +95,12 @@ class PCA * @param newDimension New dimension of the data. * @return Amount of the variance of the data retained (between 0 and 1). */ - double Apply(arma::mat& data, const size_t newDimension); + template + double Apply(MatType& data, const size_t newDimension); //! This overload is here to make sure int gets casted right to size_t. - inline double Apply(arma::mat& data, const int newDimension) + template + inline double Apply(MatType& data, const int newDimension) { return Apply(data, size_t(newDimension)); } @@ -111,7 +120,8 @@ class PCA * between 0 and 1. * @return Actual amount of variance retained (between 0 and 1). */ - double Apply(arma::mat& data, const double varRetained); + template + double Apply(MatType& data, const double varRetained); //! Get whether or not this PCA object will scale (by standard deviation) //! the data when PCA is performed. @@ -122,14 +132,15 @@ class PCA private: //! Scaling the data is when we reduce the variance of each dimension to 1. - void ScaleData(arma::mat& centeredData) + template + void ScaleData(MatType& centeredData) { if (scaleData) { // Scaling the data is when we reduce the variance of each dimension // to 1. We do this by dividing each dimension by its standard // deviation. - arma::vec stdDev = arma::stddev( + arma::Col stdDev = arma::stddev( centeredData, 0, 1 /* for each dimension */); // If there are any zeroes, make them very small. diff --git a/src/mlpack/methods/pca/pca_impl.hpp b/src/mlpack/methods/pca/pca_impl.hpp index 05db905443..a80c4486d0 100644 --- a/src/mlpack/methods/pca/pca_impl.hpp +++ b/src/mlpack/methods/pca/pca_impl.hpp @@ -37,19 +37,31 @@ PCA::PCA( * @param eigvec - PCA Loadings/Coeffs/EigenVectors */ template -void PCA::Apply(const arma::mat& data, - arma::mat& transformedData, - arma::vec& eigVal, - arma::mat& eigvec) +template +void PCA::Apply(const MatType& data, + OutMatType& transformedData, + VecType& eigVal, + OutMatType& eigvec) { + // Sanity checks on input types. + static_assert(IsBaseMatType::value, + "PCA::Apply(): transformedData must be a matrix type!"); + static_assert(IsBaseMatType::value, + "PCA::Apply(): eigVal must be a vector type!"); + static_assert(std::is_same::value, + "PCA::Apply(): data and transformedData must have the same element " + "types!"); + // Center the data into a temporary matrix. - arma::mat centeredData = data.each_col() - arma::mean(data, 1); + OutMatType centeredData = data; + centeredData.each_col() -= arma::mean(data, 1); // Scale the data if the user asked for it. ScaleData(centeredData); decomposition.Apply(data, centeredData, transformedData, eigVal, eigvec, - data.n_rows); + centeredData.n_rows); } /** @@ -60,11 +72,22 @@ void PCA::Apply(const arma::mat& data, * @param eigVal - contains eigen values in a column vector */ template -void PCA::Apply(const arma::mat& data, - arma::mat& transformedData, - arma::vec& eigVal) +template +void PCA::Apply(const MatType& data, + OutMatType& transformedData, + VecType& eigVal) { - arma::mat eigvec; + // Sanity checks on input types. + static_assert(IsBaseMatType::value, + "PCA::Apply(): transformedData must be a matrix type!"); + static_assert(IsBaseMatType::value, + "PCA::Apply(): eigVal must be a vector type!"); + static_assert(std::is_same::value, + "PCA::Apply(): data and transformedData must have the same element " + "types!"); + + OutMatType eigvec; Apply(data, transformedData, eigVal, eigvec); } @@ -75,11 +98,24 @@ void PCA::Apply(const arma::mat& data, * @param transformedData Data with PCA applied. */ template -void PCA::Apply(const arma::mat& data, - arma::mat& transformedData) +template +void PCA::Apply(const MatType& data, + OutMatType& transformedData) { - arma::mat eigvec; - arma::vec eigVal; + // Sanity checks on input types. + static_assert(IsBaseMatType::value, + "PCA::Apply(): transformedData must be a matrix type!"); + static_assert(std::is_same::value, + "PCA::Apply(): data and transformedData must have the same element " + "types!"); + + // It's possible a user didn't pass in a matrix but instead an expression, but + // we need a type that we can store. + typedef typename GetDenseColType::type BaseColType; + + OutMatType eigvec; + BaseColType eigVal; Apply(data, transformedData, eigVal, eigvec); } @@ -95,23 +131,39 @@ void PCA::Apply(const arma::mat& data, * @return Amount of the variance of the data retained (between 0 and 1). */ template -double PCA::Apply(arma::mat& data, +template +double PCA::Apply(MatType& data, const size_t newDimension) { // Parameter validation. if (newDimension == 0) - Log::Fatal << "PCA::Apply(): newDimension (" << newDimension << ") cannot " - << "be zero!" << std::endl; - if (newDimension > data.n_rows) - Log::Fatal << "PCA::Apply(): newDimension (" << newDimension << ") cannot " - << "be greater than the existing dimensionality of the data (" - << data.n_rows << ")!" << std::endl; + { + std::ostringstream oss; + oss << "PCA::Apply(): newDimension (" << newDimension << ") cannot be " + << "zero!"; + throw std::invalid_argument(oss.str()); + } - arma::mat eigvec; - arma::vec eigVal; + typedef typename GetDenseMatType::type BaseMatType; + typedef typename GetDenseColType::type BaseColType; + + BaseMatType eigvec; + BaseColType eigVal; // Center the data into a temporary matrix. - arma::mat centeredData = data.each_col() - arma::mean(data, 1); + BaseMatType centeredData = data; + centeredData.each_col() -= arma::mean(data, 1); + + // This check cannot happen until here, as `data` may not have a .n_rows + // member if it is an expression. + if (newDimension > centeredData.n_rows) + { + std::ostringstream oss; + oss << "PCA::Apply(): newDimension (" << newDimension << ") cannot " + << "be greater than the existing dimensionality of the data (" + << data.n_rows << ")!"; + throw std::invalid_argument(oss.str()); + } // Scale the data if the user ask for. ScaleData(centeredData); @@ -141,21 +193,35 @@ double PCA::Apply(arma::mat& data, * always be greater than or equal to the varRetained parameter. */ template -double PCA::Apply(arma::mat& data, +template +double PCA::Apply(MatType& data, const double varRetained) { // Parameter validation. if (varRetained < 0) - Log::Fatal << "PCA::Apply(): varRetained (" << varRetained << ") must be " - << "greater than or equal to 0." << std::endl; - if (varRetained > 1) - Log::Fatal << "PCA::Apply(): varRetained (" << varRetained << ") should be " - << "less than or equal to 1." << std::endl; + { + std::ostringstream oss; + oss << "PCA::Apply(): varRetained (" << varRetained << ") must be greater " + << "than or equal to 0."; + throw std::invalid_argument(oss.str()); + } + else if (varRetained > 1) + { + std::ostringstream oss; + oss << "PCA::Apply(): varRetained (" << varRetained << ") should be less " + << "than or equal to 1."; + throw std::invalid_argument(oss.str()); + } - arma::mat eigvec; - arma::vec eigVal; + typedef typename GetDenseMatType::type BaseMatType; + typedef typename GetDenseColType::type BaseColType; - Apply(data, data, eigVal, eigvec); + BaseMatType eigvec; + BaseColType eigVal; + BaseMatType out; + + Apply(data, out, eigVal, eigvec); + data = std::move(out); // Calculate the dimension we should keep. size_t newDimension = 0; diff --git a/src/mlpack/methods/quic_svd/quic_svd.hpp b/src/mlpack/methods/quic_svd/quic_svd.hpp index d3c8d976a5..159fe04f1e 100644 --- a/src/mlpack/methods/quic_svd/quic_svd.hpp +++ b/src/mlpack/methods/quic_svd/quic_svd.hpp @@ -51,6 +51,7 @@ namespace mlpack { * qSVD.Apply(data, u, v, sigma, epsilon, delta); * @endcode */ +template class QUIC_SVD { public: @@ -64,10 +65,10 @@ class QUIC_SVD * @param epsilon Error tolerance fraction for calculated subspace. * @param delta Cumulative probability for Monte Carlo error lower bound. */ - QUIC_SVD(const arma::mat& dataset, - arma::mat& u, - arma::mat& v, - arma::mat& sigma, + QUIC_SVD(const MatType& dataset, + MatType& u, + MatType& v, + MatType& sigma, const double epsilon = 0.03, const double delta = 0.1); @@ -93,10 +94,10 @@ class QUIC_SVD * @param epsilon Error tolerance fraction for calculated subspace. * @param delta Cumulative probability for Monte Carlo error lower bound. */ - void Apply(const arma::mat& dataset, - arma::mat& u, - arma::mat& v, - arma::mat& sigma, + void Apply(const MatType& dataset, + MatType& u, + MatType& v, + MatType& sigma, const double epsilon = 0.03, const double delta = 0.1); @@ -108,14 +109,14 @@ class QUIC_SVD * @param v Second unitary matrix. * @param sigma Diagonal matrix of singular values. */ - void ExtractSVD(const arma::mat& dataset, - arma::mat& u, - arma::mat& v, - arma::mat& sigma); + void ExtractSVD(const MatType& dataset, + MatType& u, + MatType& v, + MatType& sigma); private: //! Subspace basis of the input dataset. - arma::mat basis; + MatType basis; }; } // namespace mlpack diff --git a/src/mlpack/methods/quic_svd/quic_svd_impl.hpp b/src/mlpack/methods/quic_svd/quic_svd_impl.hpp index 4938c62892..f2a0d00124 100644 --- a/src/mlpack/methods/quic_svd/quic_svd_impl.hpp +++ b/src/mlpack/methods/quic_svd/quic_svd_impl.hpp @@ -17,39 +17,42 @@ namespace mlpack { -inline QUIC_SVD::QUIC_SVD( - const arma::mat& dataset, - arma::mat& u, - arma::mat& v, - arma::mat& sigma, +template +inline QUIC_SVD::QUIC_SVD( + const MatType& dataset, + MatType& u, + MatType& v, + MatType& sigma, const double epsilon, const double delta) { Apply(dataset, u, v, sigma, epsilon, delta); } -inline QUIC_SVD::QUIC_SVD( +template +inline QUIC_SVD::QUIC_SVD( const double /* epsilon */, const double /* delta */) { /* Nothing to do here */ } -inline void QUIC_SVD::Apply( - const arma::mat& dataset, - arma::mat& u, - arma::mat& v, - arma::mat& sigma, +template +inline void QUIC_SVD::Apply( + const MatType& dataset, + MatType& u, + MatType& v, + MatType& sigma, const double epsilon, const double delta) { // Since columns are sample in the implementation, the matrix is transposed if // necessary for maximum speedup. - CosineTree* ctree; + CosineTree* ctree; if (dataset.n_cols > dataset.n_rows) - ctree = new CosineTree(dataset, epsilon, delta); + ctree = new CosineTree(dataset, epsilon, delta); else - ctree = new CosineTree(dataset.t(), epsilon, delta); + ctree = new CosineTree(dataset.t(), epsilon, delta); // Get subspace basis by creating the cosine tree. ctree->GetFinalBasis(basis); @@ -62,24 +65,25 @@ inline void QUIC_SVD::Apply( ExtractSVD(dataset, u, v, sigma); } -inline void QUIC_SVD::ExtractSVD(const arma::mat& dataset, - arma::mat& u, - arma::mat& v, - arma::mat& sigma) +template +inline void QUIC_SVD::ExtractSVD(const MatType& dataset, + MatType& u, + MatType& v, + MatType& sigma) { // Calculate A * V_hat, necessary for further calculations. - arma::mat projectedMat; + MatType projectedMat; if (dataset.n_cols > dataset.n_rows) projectedMat = dataset.t() * basis; else projectedMat = dataset * basis; // Calculate the squared projected matrix. - arma::mat projectedMatSquared = projectedMat.t() * projectedMat; + MatType projectedMatSquared = projectedMat.t() * projectedMat; // Calculate the SVD of the above matrix. - arma::mat uBar, vBar; - arma::vec sigmaBar; + MatType uBar, vBar; + arma::Col sigmaBar; arma::svd(uBar, sigmaBar, vBar, projectedMatSquared); // Calculate the approximate SVD of the original matrix, using the SVD of the @@ -92,7 +96,7 @@ inline void QUIC_SVD::ExtractSVD(const arma::mat& dataset, // the transposed matrix is not passed. if (dataset.n_cols > dataset.n_rows) { - arma::mat tempMat = u; + MatType tempMat = u; u = v; v = tempMat; } diff --git a/src/mlpack/methods/randomized_svd/randomized_svd.hpp b/src/mlpack/methods/randomized_svd/randomized_svd.hpp index 652fe3032a..80a6782293 100644 --- a/src/mlpack/methods/randomized_svd/randomized_svd.hpp +++ b/src/mlpack/methods/randomized_svd/randomized_svd.hpp @@ -80,10 +80,11 @@ class RandomizedSVD * @param eps The eps coefficient to avoid division by zero (numerical * stability). */ - RandomizedSVD(const arma::mat& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, + template + RandomizedSVD(const InMatType& data, + MatType& u, + VecType& s, + MatType& v, const size_t iteratedPower = 0, const size_t maxIterations = 2, const size_t rank = 0, @@ -113,10 +114,11 @@ class RandomizedSVD * @param s Diagonal "Sigma" matrix of singular values. * @param rank Rank of the approximation. */ - void Apply(const arma::sp_mat& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, + template + void Apply(const arma::SpMat& data, + MatType& u, + VecType& s, + MatType& v, const size_t rank); /** @@ -129,10 +131,11 @@ class RandomizedSVD * @param s Diagonal "Sigma" matrix of singular values. * @param rank Rank of the approximation. */ - void Apply(const arma::mat& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, + template + void Apply(const InMatType& data, + MatType& u, + VecType& s, + MatType& v, const size_t rank); /** @@ -146,13 +149,16 @@ class RandomizedSVD * @param rank Rank of the approximation. * @param rowMean Centered mean value matrix. */ - template - void Apply(const MatType& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, + template + void Apply(const InMatType& data, + MatType& u, + VecType& s, + MatType& v, const size_t rank, - MatType rowMean); + const MeanType& rowMean); //! Get the size of the normalized power iterations. size_t IteratedPower() const { return iteratedPower; } diff --git a/src/mlpack/methods/randomized_svd/randomized_svd_impl.hpp b/src/mlpack/methods/randomized_svd/randomized_svd_impl.hpp index 353c68a0ea..dbc3c88eeb 100644 --- a/src/mlpack/methods/randomized_svd/randomized_svd_impl.hpp +++ b/src/mlpack/methods/randomized_svd/randomized_svd_impl.hpp @@ -17,11 +17,12 @@ namespace mlpack { +template inline RandomizedSVD::RandomizedSVD( - const arma::mat& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, + const InMatType& data, + MatType& u, + VecType& s, + MatType& v, const size_t iteratedPower, const size_t maxIterations, const size_t rank, @@ -51,42 +52,50 @@ inline RandomizedSVD::RandomizedSVD( /* Nothing to do here */ } -inline void RandomizedSVD::Apply(const arma::sp_mat& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, +template +inline void RandomizedSVD::Apply(const arma::SpMat& data, + MatType& u, + VecType& s, + MatType& v, const size_t rank) { // Center the data into a temporary matrix for sparse matrix. - arma::sp_mat rowMean = sum(data, 1) / data.n_cols; + arma::SpMat rowMean = sum(data, 1) / data.n_cols; Apply(data, u, s, v, rank, rowMean); } -inline void RandomizedSVD::Apply(const arma::mat& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, +template +inline void RandomizedSVD::Apply(const InMatType& dataIn, + MatType& u, + VecType& s, + MatType& v, const size_t rank) { // Center the data into a temporary matrix. - arma::mat rowMean = sum(data, 1) / data.n_cols + eps; + MatType data; + UnwrapAlias(data, dataIn); + + MatType rowMean = sum(data, 1) / data.n_cols + eps; Apply(data, u, s, v, rank, rowMean); } -template -inline void RandomizedSVD::Apply(const MatType& data, - arma::mat& u, - arma::vec& s, - arma::mat& v, +template +inline void RandomizedSVD::Apply(const InMatType& data, + MatType& u, + VecType& s, + MatType& v, const size_t rank, - MatType rowMean) + const MeanType& rowMean) { if (iteratedPower == 0) iteratedPower = rank + 2; - arma::mat R, Q, Qdata; + MatType R, Q, Qdata; // Apply the centered data matrix to a random matrix, obtaining Q. if (data.n_cols >= data.n_rows) @@ -97,14 +106,14 @@ inline void RandomizedSVD::Apply(const MatType& data, else { R.randn(data.n_cols, iteratedPower); - Q = (data * R) - (rowMean * (arma::ones(1, data.n_cols) * R)); + Q = (data * R) - (rowMean * (arma::ones(1, data.n_cols) * R)); } // Form a matrix Q whose columns constitute a // well-conditioned basis for the columns of the earlier Q. if (maxIterations == 0) { - arma::qr_econ(Q, v, Q); + arma::qr_econ(Q, v, Q); } else { @@ -116,7 +125,7 @@ inline void RandomizedSVD::Apply(const MatType& data, { if (data.n_cols >= data.n_rows) { - Q = (data * Q) - rowMean * (arma::ones(1, data.n_cols) * Q); + Q = (data * Q) - rowMean * (arma::ones(1, data.n_cols) * Q); arma::lu(Q, v, Q); Q = (data.t() * Q) - repmat(rowMean.t() * Q, data.n_cols, 1); } @@ -124,7 +133,7 @@ inline void RandomizedSVD::Apply(const MatType& data, { Q = (data.t() * Q) - repmat(rowMean.t() * Q, data.n_cols, 1); arma::lu(Q, v, Q); - Q = (data * Q) - (rowMean * (arma::ones(1, data.n_cols) * Q)); + Q = (data * Q) - (rowMean * (arma::ones(1, data.n_cols) * Q)); } // Computing the LU decomposition is more efficient than computing the QR @@ -146,7 +155,7 @@ inline void RandomizedSVD::Apply(const MatType& data, // applied to Q. if (data.n_cols >= data.n_rows) { - Qdata = (data * Q) - rowMean * (arma::ones(1, data.n_cols) * Q); + Qdata = (data * Q) - rowMean * (arma::ones(1, data.n_cols) * Q); arma::svd_econ(u, s, v, Qdata); v = Q * v; } diff --git a/src/mlpack/tests/cosine_tree_test.cpp b/src/mlpack/tests/cosine_tree_test.cpp index 05b807c60d..4719eaaa0a 100644 --- a/src/mlpack/tests/cosine_tree_test.cpp +++ b/src/mlpack/tests/cosine_tree_test.cpp @@ -34,7 +34,7 @@ TEST_CASE("CosineTreeNoSplit", "[CosineTreeTest]") // Make a cosine tree, with the generated dataset and the defined constants. // Note that the value of epsilon is one. - CosineTree ctree(data, epsilon, delta); + CosineTree<> ctree(data, epsilon, delta); arma::mat basis; ctree.GetFinalBasis(basis); @@ -57,17 +57,17 @@ TEST_CASE("CosineNodeCosineSplit", "[CosineTreeTest]") // Make a random dataset and the root object. arma::mat data = arma::randu(numRows, numCols); - CosineTree root(data); + CosineTree<> root(data); // Stack for depth first search of the tree. - std::vector nodeStack; + std::vector*> nodeStack; nodeStack.push_back(&root); // While stack is not empty. while (nodeStack.size()) { // Pop a node from the stack and split it. - CosineTree *currentNode, *currentLeft, *currentRight; + CosineTree<> *currentNode, *currentLeft, *currentRight; currentNode = nodeStack.back(); currentNode->CosineNodeSplit(); nodeStack.pop_back(); @@ -173,14 +173,14 @@ TEST_CASE("CosineTreeModifiedGramSchmidt", "[CosineTreeTest]") // Declare a queue and a dummy CosineTree object. CompareCosineNode comp; - CosineNodeQueue basisQueue; - CosineTree dummyTree(data, epsilon, delta); + CosineNodeQueue<> basisQueue; + CosineTree<> dummyTree(data, epsilon, delta); for (size_t i = 0; i < numCols; ++i) { // Make a new CosineNode object. - CosineTree* basisNode; - basisNode = new CosineTree(data); + CosineTree<>* basisNode; + basisNode = new CosineTree<>(data); // Use the columns of the dataset as random centroids. arma::vec centroid = data.col(i); @@ -190,8 +190,8 @@ TEST_CASE("CosineTreeModifiedGramSchmidt", "[CosineTreeTest]") dummyTree.ModifiedGramSchmidt(basisQueue, centroid, newBasisVector); // Check if the obtained vector is orthonormal to the basis vectors. - CosineNodeQueue::const_iterator j = basisQueue.cbegin(); - CosineTree* currentNode; + CosineNodeQueue<>::const_iterator j = basisQueue.cbegin(); + CosineTree<>* currentNode; for (; j != basisQueue.cend(); ++j) { @@ -210,7 +210,7 @@ TEST_CASE("CosineTreeModifiedGramSchmidt", "[CosineTreeTest]") // Deallocate memory given to the objects. for (size_t i = 0; i < numCols; ++i) { - CosineTree* currentNode; + CosineTree<>* currentNode; currentNode = basisQueue.front(); std::pop_heap(basisQueue.begin(), basisQueue.end(), comp); basisQueue.pop_back(); @@ -236,17 +236,17 @@ TEST_CASE("CopyConstructorAndOperatorCosineTreeTest", "[CosineTreeTest]") arma::mat* data = new arma::mat(numRows, numCols, arma::fill::randu); // Make a cosine tree, with the generated dataset. - CosineTree* ctree1 = new CosineTree(*data); + CosineTree<>* ctree1 = new CosineTree<>(*data); // Stacks for depth first search of the tree. - std::vector nodeStack1, nodeStack2, nodeStack3; + std::vector*> nodeStack1, nodeStack2, nodeStack3; nodeStack1.push_back(ctree1); // While stack is not empty. while (nodeStack1.size()) { // Pop a node from the stack and split it. - CosineTree *currentNode1, *currentLeft1, *currentRight1; + CosineTree<> *currentNode1, *currentLeft1, *currentRight1; currentNode1 = nodeStack1.back(); currentNode1->CosineNodeSplit(); @@ -268,8 +268,8 @@ TEST_CASE("CopyConstructorAndOperatorCosineTreeTest", "[CosineTreeTest]") } // Copy constructor and operator. - CosineTree ctree2(*ctree1); - CosineTree ctree3 = *ctree1; + CosineTree<> ctree2(*ctree1); + CosineTree<> ctree3 = *ctree1; delete ctree1; delete data; @@ -281,8 +281,8 @@ TEST_CASE("CopyConstructorAndOperatorCosineTreeTest", "[CosineTreeTest]") while (nodeStack2.size() && nodeStack3.size()) { // Pop a node from the stack and split it. - CosineTree *currentNode2, *currentLeft2, *currentRight2; - CosineTree *currentNode3, *currentLeft3, *currentRight3; + CosineTree<> *currentNode2, *currentLeft2, *currentRight2; + CosineTree<> *currentNode3, *currentLeft3, *currentRight3; currentNode2 = nodeStack2.back(); nodeStack2.pop_back(); @@ -337,17 +337,17 @@ TEST_CASE("MoveConstructorAndOperatorCosineTreeTest", "[CosineTreeTest]") arma::mat data = arma::randu(numRows, numCols); // Make a cosine tree, with the generated dataset. - CosineTree ctree1(data); + CosineTree<> ctree1(data); // Stacks for depth first search of the tree. - std::vector nodeStack1, nodeStack2, nodeStack3; + std::vector*> nodeStack1, nodeStack2, nodeStack3; nodeStack1.push_back(&ctree1); // While stack is not empty. while (nodeStack1.size()) { // Pop a node from the stack and split it. - CosineTree *currentNode1, *currentLeft1, *currentRight1; + CosineTree<> *currentNode1, *currentLeft1, *currentRight1; currentNode1 = nodeStack1.back(); currentNode1->CosineNodeSplit(); @@ -369,7 +369,7 @@ TEST_CASE("MoveConstructorAndOperatorCosineTreeTest", "[CosineTreeTest]") } // Move constructor. - CosineTree ctree2(std::move(ctree1)); + CosineTree<> ctree2(std::move(ctree1)); nodeStack2.push_back(&ctree2); @@ -377,7 +377,7 @@ TEST_CASE("MoveConstructorAndOperatorCosineTreeTest", "[CosineTreeTest]") while (nodeStack2.size()) { // Pop a node from the stack and split it. - CosineTree *currentNode2, *currentLeft2, *currentRight2; + CosineTree<> *currentNode2, *currentLeft2, *currentRight2; currentNode2 = nodeStack2.back(); nodeStack2.pop_back(); @@ -398,7 +398,7 @@ TEST_CASE("MoveConstructorAndOperatorCosineTreeTest", "[CosineTreeTest]") } // Move operator. - CosineTree ctree3 = std::move(ctree2); + CosineTree<> ctree3 = std::move(ctree2); nodeStack3.push_back(&ctree3); @@ -406,7 +406,7 @@ TEST_CASE("MoveConstructorAndOperatorCosineTreeTest", "[CosineTreeTest]") while (nodeStack3.size()) { // Pop a node from the stack and split it. - CosineTree *currentNode3, *currentLeft3, *currentRight3; + CosineTree<> *currentNode3, *currentLeft3, *currentRight3; currentNode3 = nodeStack3.back(); nodeStack3.pop_back(); diff --git a/src/mlpack/tests/pca_test.cpp b/src/mlpack/tests/pca_test.cpp index c269ac7e39..332934e662 100644 --- a/src/mlpack/tests/pca_test.cpp +++ b/src/mlpack/tests/pca_test.cpp @@ -330,3 +330,140 @@ TEST_CASE("PCAScalingTest", "[PCATest]") // The eigenvalues should sum to three. REQUIRE(accu(eigval) == Approx(3.0).epsilon(0.001)); } + +/** + * Test PCA on a subview of a matrix with different decomposition strategies. + */ +TEMPLATE_TEST_CASE("PCASubviewTest", "[PCATest]", ExactSVDPolicy, + RandomizedSVDPCAPolicy, RandomizedBlockKrylovSVDPolicy, QUICSVDPolicy) +{ + typedef TestType DecompositionPolicy; + + // Generate an artifical dataset in 10 dimensions. + arma::mat data(3, 5000); + + arma::vec mean("1.0 3.0 -12.0"); + arma::mat cov("1.0 0.9 0.0;" + "0.9 1.0 0.0;" + "0.0 0.0 12.0"); + GaussianDistribution g(mean, cov); + + for (size_t i = 0; i < 5000; ++i) + data.col(i) = g.Random(); + + // Compute PCA on the first 2000 points. + arma::mat transData1, transData2, transData3, eigvec; + arma::vec eigval1, eigval2; + + PCA p; + p.Apply(data.cols(0, 1999), transData1); + p.Apply(data.cols(0, 1999), transData2, eigval1); + p.Apply(data.cols(0, 1999), transData3, eigval2, eigvec); + + // Only check for deterministic policies. + if (std::is_same::value) + { + arma::mat trueTransData, trueEigvec; + arma::vec trueEigval; + + arma::mat dataSub = data.cols(0, 1999); + p.Apply(dataSub, trueTransData, trueEigval, trueEigvec); + + REQUIRE(arma::approx_equal(transData1, trueTransData, "both", 1e-5, 1e-5)); + REQUIRE(arma::approx_equal(transData2, trueTransData, "both", 1e-5, 1e-5)); + REQUIRE(arma::approx_equal(transData3, trueTransData, "both", 1e-5, 1e-5)); + REQUIRE(arma::approx_equal(eigval1, trueEigval, "both", 1e-5, 1e-5)); + REQUIRE(arma::approx_equal(eigval2, trueEigval, "both", 1e-5, 1e-5)); + REQUIRE(arma::approx_equal(eigvec, trueEigvec, "both", 1e-5, 1e-5)); + } +} + +/** + * Test PCA on an input expression. + */ +TEMPLATE_TEST_CASE("PCAExpressionTest", "[PCATest]", ExactSVDPolicy, + RandomizedSVDPCAPolicy, RandomizedBlockKrylovSVDPolicy, QUICSVDPolicy) +{ + typedef TestType DecompositionPolicy; + + // Generate an artifical dataset in 10 dimensions. + arma::mat data(3, 5000); + + arma::vec mean("1.0 3.0 -12.0"); + arma::mat cov("1.0 0.9 0.0;" + "0.9 1.0 0.0;" + "0.0 0.0 12.0"); + GaussianDistribution g(mean, cov); + + for (size_t i = 0; i < 5000; ++i) + data.col(i) = g.Random(); + + // Compute PCA on an expression involving the input matrix. + arma::mat transData1, transData2, transData3, eigvec; + arma::vec eigval1, eigval2; + + PCA p; + p.Apply(2 * data + 1, transData1); + p.Apply(2 * data + 1, transData2, eigval1); + p.Apply(2 * data + 1, transData3, eigval2, eigvec); + + // Only check for deterministic policies. + if (std::is_same::value) + { + arma::mat trueTransData, trueEigvec; + arma::vec trueEigval; + + arma::mat dataSub = 2 * data + 1; + p.Apply(dataSub, trueTransData, trueEigval, trueEigvec); + + REQUIRE(arma::approx_equal(transData1, trueTransData, "both", 1e-5, 1e-5)); + REQUIRE(arma::approx_equal(transData2, trueTransData, "both", 1e-5, 1e-5)); + REQUIRE(arma::approx_equal(transData3, trueTransData, "both", 1e-5, 1e-5)); + REQUIRE(arma::approx_equal(eigval1, trueEigval, "both", 1e-5, 1e-5)); + REQUIRE(arma::approx_equal(eigval2, trueEigval, "both", 1e-5, 1e-5)); + REQUIRE(arma::approx_equal(eigvec, trueEigvec, "both", 1e-5, 1e-5)); + } +} + +/** + * Test PCA on 32-bit data. + */ +TEMPLATE_TEST_CASE("PCAFloatTest", "[PCATest]", ExactSVDPolicy, + RandomizedSVDPCAPolicy, RandomizedBlockKrylovSVDPolicy, QUICSVDPolicy) +{ + typedef TestType DecompositionPolicy; + + // Generate an artifical dataset in 10 dimensions. + arma::fmat data(3, 5000); + + arma::vec mean("1.0 3.0 -12.0"); + arma::mat cov("1.0 0.9 0.0;" + "0.9 1.0 0.0;" + "0.0 0.0 12.0"); + GaussianDistribution g(mean, cov); + + for (size_t i = 0; i < 5000; ++i) + data.col(i) = arma::conv_to::from(g.Random()); + + // Compute PCA on the floating-point data. + arma::fmat coeff, coeff1, score, score1; + arma::fvec eigVal, eigVal1; + + PCA p; + p.Apply(data, score1, eigVal1, coeff1); + + princomp(coeff, score, eigVal, trans(data)); + + // Verify the PCA results based on the eigenvalues. We don't check for + // QUIC-SVD, since that method has a lot of noise. + if (!std::is_same::value) + { + for (size_t i = 0; i < eigVal.n_elem; ++i) + { + if (eigVal[i] == 0.0) + REQUIRE(eigVal1[i] == Approx(0.0).margin(1e-5)); + else + REQUIRE(eigVal[i] == Approx(eigVal1[i]).epsilon(1e-3)); + } + } +} diff --git a/src/mlpack/tests/policy_gradient_test.cpp b/src/mlpack/tests/policy_gradient_test.cpp index f07a710c79..7b92f78a9e 100644 --- a/src/mlpack/tests/policy_gradient_test.cpp +++ b/src/mlpack/tests/policy_gradient_test.cpp @@ -192,7 +192,7 @@ TEST_CASE("OUNoiseTest", "[PolicyGradientTest]") // Test the reset function. ouNoise.reset(); arma::colvec state = ouNoise.sample(); - REQUIRE(state.n_elem == size); + REQUIRE((int) state.n_elem == size); // Verify that the sample is not equal to the reset state. arma::colvec sample = ouNoise.sample(); @@ -213,7 +213,7 @@ TEST_CASE("GaussianNoiseTest", "[PolicyGradientTest]") // Test the sample function. arma::colvec noise = gaussianNoise.sample(); - REQUIRE(noise.n_elem == size); + REQUIRE((int) noise.n_elem == size); // Verify that the noise vector has values drawn from a // Gaussian distribution with the specified mean and standard deviation. diff --git a/src/mlpack/tests/quic_svd_test.cpp b/src/mlpack/tests/quic_svd_test.cpp index b77fca3750..772172c7b1 100644 --- a/src/mlpack/tests/quic_svd_test.cpp +++ b/src/mlpack/tests/quic_svd_test.cpp @@ -34,7 +34,7 @@ TEST_CASE("QUICSVDReconstructionError", "[QUICSVDTest]") { // Obtain the SVD using default parameters. arma::mat u, v, sigma; - QUIC_SVD quicsvd(dataset, u, v, sigma); + QUIC_SVD<> quicsvd(dataset, u, v, sigma); // Reconstruct the matrix using the SVD. arma::mat reconstruct; @@ -71,7 +71,7 @@ TEST_CASE("QUICSVDSingularValueError", "[QUICSVDTest]") // Obtain the SVD using default parameters. arma::svd_econ(U1, s1, V1, data); - QUIC_SVD quicsvd(data, U1, V1, s2); + QUIC_SVD<> quicsvd(data, U1, V1, s2); s3 = arma::diagvec(s2); s1 = s1.subvec(0, s3.n_elem - 1); @@ -87,5 +87,5 @@ TEST_CASE("QUICSVDSameDimensionTest", "[QUICSVDTest]") // Obtain the SVD using default parameters. arma::mat u, v, sigma; - QUIC_SVD quicsvd(dataset, u, v, sigma); + QUIC_SVD<> quicsvd(dataset, u, v, sigma); } From b3fd17e65fb25d3644b7b452a46df39782cf1678 Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Tue, 2 Apr 2024 14:54:21 -0400 Subject: [PATCH 03/11] Add documentation for PCA. --- doc/index.md | 2 +- doc/sidebar.html | 17 ++- doc/user/methods/pca.md | 284 ++++++++++++++++++++++++++++++++++++++++ 3 files changed, 299 insertions(+), 4 deletions(-) create mode 100644 doc/user/methods/pca.md diff --git a/doc/index.md b/doc/index.md index e8d64e4d72..e799782e41 100644 --- a/doc/index.md +++ b/doc/index.md @@ -120,7 +120,7 @@ Prepare data for machine learning algorithms. Transform data from one space to another. - + * [`PCA`](pca.md): principal components analysis ### Modeling utilities diff --git a/doc/sidebar.html b/doc/sidebar.html index da95a2d834..9890759b00 100644 --- a/doc/sidebar.html +++ b/doc/sidebar.html @@ -135,9 +135,20 @@ when the sidebar is built for each page.
  • - - Transformations - +
    + + + Transformations + + + +
    • diff --git a/doc/user/methods/pca.md b/doc/user/methods/pca.md new file mode 100644 index 0000000000..e58605803d --- /dev/null +++ b/doc/user/methods/pca.md @@ -0,0 +1,284 @@ +## `PCA` + +The `PCA` class implements principal components analysis (PCA), a standard +machine learning data preparation technique. PCA can be used to reduce the +number of dimensions in a dataset, or to preseve a certain percentage of the +variance of a dataset. + +By default, `PCA` uses the singular value decomposition (SVD), but supports the +use of other more efficient decompositions, including approximate +decompositions. + +#### Simple usage example: + +```c++ +// Use PCA to reduce the number of dimensions to 5 on uniform random data. + +// This dataset is uniform random in 10 dimensions. +// Replace with a data::Load() call or similar for a real application. +arma::mat dataset(10, 1000, arma::fill::randu); // 1000 points. + +mlpack::PCA pca; // Step 1: create PCA object. +pca.Apply(data, 5); // Step 2: reduce data dimension to 5. + +// Print some information about the modified dataset. +std::cout << "The transformed data matrix has size " << data.n_rows /* 5 */ + << " x " << data.n_cols << "." << std::endl; +``` +

    More examples...

    + +#### Quick links: + + * [Constructors](#constructors): create `PCA` objects. + * [`Apply()`](#applying-transformation): apply PCA transformation to data. + * [Examples](#simple-examples) of simple usage and links to detailed example + projects. + * [Template parameters](#advanced-functionality-different-decomposition-strategies) + for using different decomposition strategies. + +#### See also: + + * [`RADICAL`](radical.md): independent components analysis + * [mlpack preprocessing utilities](index.md#preprocessing-utilities) + * [mlpack transformations](index.md#transformations) + * [Principal component analysis on Wikipedia](https://en.wikipedia.org/wiki/Principal_component_analysis) + +### Constructors + + * `pca = PCA(scaleData=false)` + - Construct a `PCA` object. + - If `scaleData` is `true`, then all dimensions will have variance scaled to + 1 before applying PCA. + - The `scaleData` parameter can be inspected with `pca.ScaleData()`, and also + set; `pca.ScaleData() = true` will enable data variance scaling. + +--- + + * `pca = PCA(scaleData, decompositionPolicy)` + - Construct a `PCA` object with a custom decomposition policy. + - See the documentation for using + [different decomposition strategies](#advanced-functionality-different-decomposition-strategies). + +### Applying Transformations + + * `pca.Apply(data, transformedData)` + * `pca.Apply(data, transformedData, eigVal)` + * `pca.Apply(data, transformedData, eigVal, eigVec)` + - Transform the + [column-major matrix](../matrices.md#representing-data-in-mlpack) `data` + using PCA, storing the result in `transformedData`. + - `data` and `transformedData` should be dense floating-point matrices (e.g., + `arma::mat`, `arma::fmat`, etc.). + - The size of `transformedData` will be the same as the size of `data`. + - Dimensions in `transformedData` will be ordered decreasing in variance; + that is, the first row of `transformedData` will correspond to the + dimension with maximum variance. + - Optionally, eigenvalues and eigenvectors of the covariance matrix can be + returned: + * If specified, `eigVal` should be a dense floating-point vector (e.g. + `arma::vec`, `arma::fvec`, etc.) and will be filled with the eigenvalues + of `transformedData`. + * If specified, `eigvec` should be a dense floating-point matrix (e.g. + `arma::mat`, `arma::fmat`, etc.) and will be filled with the eigenvectors + of `transformedData`. + +--- + + * `double varRetained = pca.Apply(data, newDimension)` + - Use PCA to reduce the number of dimensions in the + [column-major matrix](../matrices.md#representing-data-in-mlpack) `data` + to `newDimension`. + - `data` should be a dense floating-point matrix (e.g. `arma::mat`, + `arma::fmat`, etc.). + - `data` will be modified in-place, and will have `newDimension` rows after + the transformation. + - Returns a `double` indicating the percentage of variance retained (between + `0.0` and `1.0`). + +--- + + * `double varRetained = pca.Apply(data, varianceToKeep)` + - Use PCA to retain the dimensions of the + [column-major matrix](../matrices.md#representing-data-in-mlpack) `data` + that capture a factor of `varianceToKeep` of the data variance. + - `data` should be a dense floating-point matrix (e.g. `arma::mat`, + `arma::fmat`, etc.). + - `varianceToKeep` should be a floating-point value between `0.0` and `1.0`. + If `1.0`, all of the data variance is retained, and this is equivalent to + the first version of `Apply()` (above). + - `data` will be modified in-place. + - Returns a `double` indicating the percentage of variance actually retained + (between `0.0` and `1.0`). + +--- + +### Simple Examples + +See also the [simple usage example](#simple-usage-example) for a trivial usage +of the `PCA` class. + +--- + +Apply PCA to a dataset, keeping dimensions that capture 90% of the data +variance. + +```c++ +// See https://datasets.mlpack.org/satellite.train.csv. +arma::mat data; +mlpack::data::Load("satellite.train.csv", data, true); +const size_t origDim = data.n_rows; + +mlpack::PCA pca; + +// Keep 90% of the data variance. +pca.Apply(data, 0.9); + +std::cout << "PCA kept " << data.n_rows << " of " << origDim << " dimensions " + << "to capture 90\% of the data variance." << std::endl; +``` + +--- + +Apply PCA to a 32-bit floating point dataset with dimension scaling, keeping all +dimensions, and printing the 5 largest eigenvalues of the covariance matrix of +the transformed data. + +``` +// See https://datasets.mlpack.org/iris.csv +arma::fmat data; +mlpack::data::Load("iris.csv", data, true); + +mlpack::PCA pca(true /* scale data when transforming */); + +arma::fvec eigval; +arma::fmat transformedData; + +pca.Apply(data, transformedData, eigval); + +std::cout << "First point, before PCA: " << data.col(0).t(); +std::cout << "First point, after PCA: " << transformedData.col(0).t(); +std::cout << std::endl; + +// Now print the top 5 eigenvalues. +for (size_t i = 0; i < 5; ++i) + std::cout << "Eigenvalue " << i << ": " << eigval[i] << "." << std::endl; +``` + +--- + +### Advanced Functionality: Different Decomposition Strategies + +By default, `PCA` uses the standard singular value decomposition (SVD) to +transform data. However, for very large datasets, it may be faster to use +alternative strategies, some of which may be approximate. The `PCA` class has +one template parameter that allows different decomposition strategies to be +used. The full signature of the class is: + +``` +PCA +``` + +`DecompositionPolicy` specifies the strategy to be used to compute the singular +values and vectors of a data matrix. + +Several decomposition policies are already implemented and ready for drop-in +usage: + + * `ExactSVDPolicy` _(default)_: use Armadillo's `svd()` and `svd_econ()` + functions to compute the SVD + * `RandomizedSVDPCAPolicy`: use the randomized SVD algorithm to compute the SVD + + * `RandomizedBlockKrylovSVDPolicy`: use the randomized Block Krylov SVD + algorithm to compute the SVD + * `QUICSVDPolicy`: use the tree-based `QUIC-SVD` algorithm to compute the SVD + + +The simple example program below uses all four decomposition types on the same +data, timing how long each decomposition takes. + +```c++ +arma::mat data; +// See https://datasets.mlpack.org/mnist.train.csv. +mlpack::data::Load("mnist.train.csv", true); + +arma::mat output1, output2, output3, output4; + +mlpack::PCA pca1; +mlpack::PCA pca2; +mlpack::PCA pca3; +mlpack::PCA pca4; + +// Compute decompositions on all four, timing each one. +arma::wall_clock c; + +c.tic(); +pca1.Apply(data, output1); +const double pca1Time = c.toc(); + +c.tic(); +pca2.Apply(data, output2); +const double pca2Time = c.toc(); + +c.tic(); +pca3.Apply(data, output3); +const double pca3Time = c.toc(); + +c.tic(); +pca4.Apply(data, output4); +const double pca4Time = c.toc(); + +std::cout << "PCA computation times:" << std::endl; +std::cout << " - ExactSVDPolicy: " << pca1Time << "s." + << std::endl; +std::cout << " - RandomizedSVDPCAPolicy: " << pca2Time << "s." + << std::endl; +std::cout << " - RandomizedBlockKrylovSVDPolicy: " << pca3Time << "s." + << std::endl; +std::cout << " - QUICSVDPolicy: " << pca4Time << "s." + << std::endl; +``` + +--- + +#### Custom decomposition policies + +Instead of using the predefined classes above, it is also possible to implement +fully custom functionality via a new decomposition policy. Any new +decomposition policy must implement one method: + +```c++ +class CustomDecompositionPolicy +{ + public: + // Given input data `data` and `centeredData`, compute the singular value + // decomposition of the data, and then project the data onto the first `rank` + // singular vectors. + // + // * `data` is the input matrix. It is not guaranteed to be centered or + // scaled. + // * `centeredData` is the centered (and possibly scaled) version of the + // input matrix (e.g. the mean of each dimension is 0). + // * `transformedData` should be overwritten with the centered data's + // projection onto the singular vectors. + // * `svals` and `svecs` should be filled with the singular values and + // vectors of the centered data. + // * `rank` specifies the number of singular values/vectors to keep, and the + // dimension of `transformedData` should be equivalent to `rank`. `rank` + // will be at most equal to `data.n_rows`. + // + // * `InMatType` is a dense floating-point matrix type, but may be a subview + // or expression. + // * `MatType` is the type of matrix used to represent data, and will be a + // dense floating-point matrix type (e.g. `arma::mat`, `arma::fmat`, + // etc.). + // * `VecType` is the corresponding vector type to `MatType` (e.g., a + // `MatType` of `arma::mat` would mean a `VecType` of `arma::vec`, etc.). + template + static void Apply(const InMatType& data, + const MatType& centeredData, + MatType& transformedData, + VecType& svals, + MatType& svecs, + const size_t rank); +}; +``` From 294d80da3573ffec438959b63390573fd9c511fe Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Wed, 3 Apr 2024 13:21:16 -0400 Subject: [PATCH 04/11] Revisions to PCA documentation. --- doc/user/methods/pca.md | 93 +++++++++++++++++++++++++++++------------ 1 file changed, 67 insertions(+), 26 deletions(-) diff --git a/doc/user/methods/pca.md b/doc/user/methods/pca.md index e58605803d..b01a8404f8 100644 --- a/doc/user/methods/pca.md +++ b/doc/user/methods/pca.md @@ -2,12 +2,12 @@ The `PCA` class implements principal components analysis (PCA), a standard machine learning data preparation technique. PCA can be used to reduce the -number of dimensions in a dataset, or to preseve a certain percentage of the +number of dimensions in a dataset, or to preserve a certain percentage of the variance of a dataset. -By default, `PCA` uses the singular value decomposition (SVD), but supports the -use of other more efficient decompositions, including approximate -decompositions. +By default, `PCA` uses the full exact singular value decomposition (SVD), but +supports the use of other more efficient decompositions, including approximate +singular value decompositions. #### Simple usage example: @@ -18,12 +18,12 @@ decompositions. // Replace with a data::Load() call or similar for a real application. arma::mat dataset(10, 1000, arma::fill::randu); // 1000 points. -mlpack::PCA pca; // Step 1: create PCA object. -pca.Apply(data, 5); // Step 2: reduce data dimension to 5. +mlpack::PCA pca; // Step 1: create PCA object. +pca.Apply(dataset, 5); // Step 2: reduce data dimension to 5. // Print some information about the modified dataset. -std::cout << "The transformed data matrix has size " << data.n_rows /* 5 */ - << " x " << data.n_cols << "." << std::endl; +std::cout << "The transformed data matrix has size " << dataset.n_rows /* 5 */ + << " x " << dataset.n_cols << "." << std::endl; ```

    More examples...

    @@ -67,8 +67,10 @@ std::cout << "The transformed data matrix has size " << data.n_rows /* 5 */ - Transform the [column-major matrix](../matrices.md#representing-data-in-mlpack) `data` using PCA, storing the result in `transformedData`. - - `data` and `transformedData` should be dense floating-point matrices (e.g., - `arma::mat`, `arma::fmat`, etc.). + - `data` should be a floating-point matrix (e.g. `arma::mat`, `arma::fmat`, + `arma::sp_mat`, etc.) or an expression that evaluates to one. + - `transformedData` should be a dense floating-point matrix (e.g., + `arma::mat`, `arma::sp_mat`). - The size of `transformedData` will be the same as the size of `data`. - Dimensions in `transformedData` will be ordered decreasing in variance; that is, the first row of `transformedData` will correspond to the @@ -84,32 +86,46 @@ std::cout << "The transformed data matrix has size " << data.n_rows /* 5 */ --- - * `double varRetained = pca.Apply(data, newDimension)` + * `double varRetained = pca.Apply(data, transformedData, newDimension)` - Use PCA to reduce the number of dimensions in the [column-major matrix](../matrices.md#representing-data-in-mlpack) `data` - to `newDimension`. - - `data` should be a dense floating-point matrix (e.g. `arma::mat`, - `arma::fmat`, etc.). - - `data` will be modified in-place, and will have `newDimension` rows after - the transformation. + to `newDimension`, storing the result in `transformedData`. + - `data` should be a floating-point matrix (e.g. `arma::mat`, + `arma::fmat`, `arma::sp_mat`, etc.) or an expression that evaluates to + one. + - `transformedData` should be a dense floating-point matrix with the same + element type as `data` (e.g. `arma::mat`, `arma::fmat`). + - `transformedData` will have `newDimension` rows after the transformation. - Returns a `double` indicating the percentage of variance retained (between `0.0` and `1.0`). --- - * `double varRetained = pca.Apply(data, varianceToKeep)` + * `double varRetained = pca.Apply(data, transformedData, varianceToKeep)` - Use PCA to retain the dimensions of the [column-major matrix](../matrices.md#representing-data-in-mlpack) `data` that capture a factor of `varianceToKeep` of the data variance. - - `data` should be a dense floating-point matrix (e.g. `arma::mat`, - `arma::fmat`, etc.). + - `data` should be a floating-point matrix (e.g. `arma::mat`, `arma::fmat`, + `arma::sp_mat`, etc.) or an expression that evaluates to one. + - `transformedData` should be a dense floating-point matrix with the same + element type as `data` (e.g. `arma::mat`, `arma::fmat`). + - `transformedData` will have `newDimension` rows after the transformation. - `varianceToKeep` should be a floating-point value between `0.0` and `1.0`. If `1.0`, all of the data variance is retained, and this is equivalent to the first version of `Apply()` (above). - - `data` will be modified in-place. - Returns a `double` indicating the percentage of variance actually retained (between `0.0` and `1.0`). +--- + + * `double varRetained = pca.Apply(data, newDimension)` + * `double varRetained = pca.Apply(data, varianceToKeep)` + - In-place versions of the two `Apply()` functions above. + - Equivalent to `pca.Apply(data, data, newDimension)` or + `pca.Apply(data, data, varianceToKeep)`. + - `data` should be a dense floating-point matrix (e.g. `arma::mat`, + `arma::fmat`, etc.). + --- ### Simple Examples @@ -143,8 +159,8 @@ Apply PCA to a 32-bit floating point dataset with dimension scaling, keeping all dimensions, and printing the 5 largest eigenvalues of the covariance matrix of the transformed data. -``` -// See https://datasets.mlpack.org/iris.csv +```c++ +// See https://datasets.mlpack.org/iris.csv. arma::fmat data; mlpack::data::Load("iris.csv", data, true); @@ -166,9 +182,33 @@ for (size_t i = 0; i < 5; ++i) --- +Apply PCA to a random sparse dataset, to reduce the dimensionality to a +20-dimensional dense dataset. + +```c++ +arma::sp_mat data; +// This dataset has 10k points in 1k dimensions, with 1% density. +data.sprandn(1000, 10000, 0.01); + +mlpack::PCA pca(true /* scale data when transforming */); + +arma::mat transformedData; +const double varianceRetained = pca.Apply(data, transformedData, 20); + +std::cout << "First point, before PCA: " << data.col(0).t(); +std::cout << "First point, after PCA: " << transformedData.col(0).t(); + +// Note that for random uniform data, this won't capture very much of the +// variance! It would be much more for a real, structured dataset. +std::cout << "50 dimensions captured " << (100.0 * varianceRetained) << "\% of " + << "the data variance." << std::endl; +``` + +--- + ### Advanced Functionality: Different Decomposition Strategies -By default, `PCA` uses the standard singular value decomposition (SVD) to +By default, `PCA` uses the full exact singular value decomposition (SVD) to transform data. However, for very large datasets, it may be faster to use alternative strategies, some of which may be approximate. The `PCA` class has one template parameter that allows different decomposition strategies to be @@ -194,12 +234,12 @@ usage: The simple example program below uses all four decomposition types on the same -data, timing how long each decomposition takes. +MNIST data, timing how long each decomposition takes. ```c++ arma::mat data; // See https://datasets.mlpack.org/mnist.train.csv. -mlpack::data::Load("mnist.train.csv", true); +mlpack::data::Load("mnist.train.csv", data, true); arma::mat output1, output2, output3, output4; @@ -227,7 +267,8 @@ c.tic(); pca4.Apply(data, output4); const double pca4Time = c.toc(); -std::cout << "PCA computation times:" << std::endl; +std::cout << "PCA computation times for " << data.n_rows << " x " << data.n_cols + << " data:" << std::endl; std::cout << " - ExactSVDPolicy: " << pca1Time << "s." << std::endl; std::cout << " - RandomizedSVDPCAPolicy: " << pca2Time << "s." From 2c390b0f689175d5f65335b80add3d829c2085c2 Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Wed, 3 Apr 2024 13:21:27 -0400 Subject: [PATCH 05/11] Limit width of sidebar. --- doc/css/gfm-mod.css | 1 + 1 file changed, 1 insertion(+) diff --git a/doc/css/gfm-mod.css b/doc/css/gfm-mod.css index 5cd333d7a9..29518d3836 100644 --- a/doc/css/gfm-mod.css +++ b/doc/css/gfm-mod.css @@ -1062,6 +1062,7 @@ div#sidebar { top: 5px; min-width: 200px; font-size: 90%; + max-width: 15.1515%; } div#sidebar ul { From e4f33e20c09e05a548ce21bc9e2d51b694b40127 Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Wed, 3 Apr 2024 13:21:50 -0400 Subject: [PATCH 06/11] Add additional overloads for PCA::Apply() and tests. --- .../randomized_block_krylov_svd_impl.hpp | 7 ++- src/mlpack/methods/pca/pca.hpp | 47 +++++++++++++++++++ src/mlpack/methods/pca/pca_impl.hpp | 37 +++++++++++---- src/mlpack/tests/pca_test.cpp | 40 ++++++++++++++++ 4 files changed, 119 insertions(+), 12 deletions(-) diff --git a/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp b/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp index a042546252..7a8000ebc0 100644 --- a/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp +++ b/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp @@ -58,7 +58,10 @@ inline void RandomizedBlockKrylovSVD::Apply(const InMatType& data, if (blockSize == 0) { - blockSize = rank + 10; + // The block size cannot be greater than the number of points in the + // dataset or the dimensionality of the dataset. + blockSize = std::min((size_t) data.n_rows, + std::min((size_t) data.n_cols, rank + 10)); } // Random block initialization. @@ -68,7 +71,7 @@ inline void RandomizedBlockKrylovSVD::Apply(const InMatType& data, MatType K(data.n_rows, blockSize * (maxIterations + 1)); // Create a working matrix using data from writable auxiliary memory - // (K matrix). Doing so avoids an uncessary copy in upcoming step. + // (K matrix). Doing so avoids an unnecessary copy in upcoming step. MakeAlias(block, K.memptr(), data.n_rows, blockSize, false); arma::qr_econ(block, R, data * G); diff --git a/src/mlpack/methods/pca/pca.hpp b/src/mlpack/methods/pca/pca.hpp index 612fdebbbf..bf92350d45 100644 --- a/src/mlpack/methods/pca/pca.hpp +++ b/src/mlpack/methods/pca/pca.hpp @@ -105,6 +105,33 @@ class PCA return Apply(data, size_t(newDimension)); } + /** + * Use PCA for dimensionality reduction on the given dataset. This will save + * the newDimension largest principal components of the data and remove the + * rest, storing the result in `transformedData`. The return value is the + * amount of variance of the data that is retained; this is a value between 0 + * and 1. For instance, a value of 0.9 indicates that 90% of the variance + * present in the data was retained. + * + * @param data Data matrix. + * @param transformedData Output matrix to store transformed data in. + * @param newDimension New dimension of the data. + * @return Amount of the variance of the data retained (between 0 and 1). + */ + template + double Apply(const MatType& data, + OutMatType& transformedData, + const size_t newDimension); + + //! This overload is here to make sure int gets casted right to size_t. + template + inline double Apply(const MatType& data, + OutMatType& transformedData, + const int newDimension) + { + return Apply(data, transformedData, size_t(newDimension)); + } + /** * Use PCA for dimensionality reduction on the given dataset. This will save * as many dimensions as necessary to retain at least the given amount of @@ -123,6 +150,26 @@ class PCA template double Apply(MatType& data, const double varRetained); + /** + * Use PCA for dimensionality reduction on the given dataset. This will save + * as many dimensions as necessary to retain at least the given amount of + * variance (specified by parameter varRetained). The amount should be + * between 0 and 1; if the amount is 0, then only 1 dimension will be + * retained. If the amount is 1, then all dimensions will be retained. + * + * The method returns the actual amount of variance retained, which will + * always be greater than or equal to the varRetained parameter. + * + * @param data Data matrix. + * @param varRetained Lower bound on amount of variance to retain; should be + * between 0 and 1. + * @return Actual amount of variance retained (between 0 and 1). + */ + template + double Apply(const MatType& data, + OutMatType& transformedData, + const double varRetained); + //! Get whether or not this PCA object will scale (by standard deviation) //! the data when PCA is performed. bool ScaleData() const { return scaleData; } diff --git a/src/mlpack/methods/pca/pca_impl.hpp b/src/mlpack/methods/pca/pca_impl.hpp index a80c4486d0..865c968300 100644 --- a/src/mlpack/methods/pca/pca_impl.hpp +++ b/src/mlpack/methods/pca/pca_impl.hpp @@ -12,7 +12,6 @@ * 3-clause BSD license along with mlpack. If not, see * http://www.opensource.org/licenses/BSD-3-Clause for more information. */ - #ifndef MLPACK_METHODS_PCA_PCA_IMPL_HPP #define MLPACK_METHODS_PCA_PCA_IMPL_HPP @@ -54,8 +53,8 @@ void PCA::Apply(const MatType& data, "types!"); // Center the data into a temporary matrix. - OutMatType centeredData = data; - centeredData.each_col() -= arma::mean(data, 1); + OutMatType centeredData = arma::conv_to::from(data); + centeredData.each_col() -= arma::mean(centeredData, 1); // Scale the data if the user asked for it. ScaleData(centeredData); @@ -134,6 +133,15 @@ template template double PCA::Apply(MatType& data, const size_t newDimension) +{ + return Apply(data, data, newDimension); +} + +template +template +double PCA::Apply(const MatType& data, + OutMatType& transformedData, + const size_t newDimension) { // Parameter validation. if (newDimension == 0) @@ -151,8 +159,8 @@ double PCA::Apply(MatType& data, BaseColType eigVal; // Center the data into a temporary matrix. - BaseMatType centeredData = data; - centeredData.each_col() -= arma::mean(data, 1); + BaseMatType centeredData = arma::conv_to::from(data); + centeredData.each_col() -= arma::mean(centeredData, 1); // This check cannot happen until here, as `data` may not have a .n_rows // member if it is an expression. @@ -168,11 +176,12 @@ double PCA::Apply(MatType& data, // Scale the data if the user ask for. ScaleData(centeredData); - decomposition.Apply(data, centeredData, data, eigVal, eigvec, newDimension); + decomposition.Apply(data, centeredData, transformedData, eigVal, eigvec, + newDimension); if (newDimension < eigvec.n_rows) // Drop unnecessary rows. - data.shed_rows(newDimension, data.n_rows - 1); + transformedData.shed_rows(newDimension, data.n_rows - 1); // The svd method returns only non-zero eigenvalues so we have to calculate // the right dimension before calculating the amount of variance retained. @@ -196,6 +205,15 @@ template template double PCA::Apply(MatType& data, const double varRetained) +{ + return Apply(data, data, varRetained); +} + +template +template +double PCA::Apply(const MatType& data, + OutMatType& transformedData, + const double varRetained) { // Parameter validation. if (varRetained < 0) @@ -220,8 +238,7 @@ double PCA::Apply(MatType& data, BaseColType eigVal; BaseMatType out; - Apply(data, out, eigVal, eigvec); - data = std::move(out); + Apply(data, transformedData, eigVal, eigvec); // Calculate the dimension we should keep. size_t newDimension = 0; @@ -235,7 +252,7 @@ double PCA::Apply(MatType& data, // varSum is the actual variance we will retain. if (newDimension < eigVal.n_elem) - data.shed_rows(newDimension, data.n_rows - 1); + transformedData.shed_rows(newDimension, transformedData.n_rows - 1); return varSum; } diff --git a/src/mlpack/tests/pca_test.cpp b/src/mlpack/tests/pca_test.cpp index 332934e662..c37f3c3881 100644 --- a/src/mlpack/tests/pca_test.cpp +++ b/src/mlpack/tests/pca_test.cpp @@ -467,3 +467,43 @@ TEMPLATE_TEST_CASE("PCAFloatTest", "[PCATest]", ExactSVDPolicy, } } } + +/** + * Test that we can convert sparse input matrices to dense output matrices. + * We check the general shape of the output, but not the exact details---those + * are handled in other tests. + */ +TEMPLATE_TEST_CASE("PCASparseToDenseTest", "[PCATest]", float, double) +{ + typedef arma::Mat MatType; + typedef arma::SpMat SpMatType; + + SpMatType dataset; + dataset.sprandu(1000, 50000, 0.01); + MatType transformedDataset1, transformedDataset2; + + PCA<> p; + const double varRetained1 = p.Apply(dataset, transformedDataset1, 5); + const double varRetained2 = p.Apply(dataset, transformedDataset2, 0.6); + + REQUIRE(transformedDataset1.n_cols == dataset.n_cols); + REQUIRE(transformedDataset2.n_cols == dataset.n_cols); + REQUIRE(varRetained1 >= 0.0); + REQUIRE(varRetained1 <= 1.0); + REQUIRE(varRetained2 >= 0.0); + REQUIRE(varRetained2 <= 1.0); + REQUIRE(transformedDataset1.n_rows == 5); + REQUIRE(transformedDataset2.n_rows <= dataset.n_rows); + + // Ensure we get basically the same as if we had done it to a dense matrix. + MatType denseData1(dataset); + MatType denseData2(dataset); + + const double varRetained3 = p.Apply(denseData1, 5); + const double varRetained4 = p.Apply(denseData2, 0.6); + + REQUIRE(varRetained1 == Approx(varRetained3)); + REQUIRE(varRetained2 == Approx(varRetained4)); + + REQUIRE(denseData2.n_rows == transformedDataset2.n_rows); +} From 276dc85aac46d9c3041e89e26a59ede0faffca78 Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Wed, 3 Apr 2024 13:57:34 -0400 Subject: [PATCH 07/11] Fix links in PCA documentation. --- doc/index.md | 2 +- doc/user/methods/pca.md | 9 +++++---- 2 files changed, 6 insertions(+), 5 deletions(-) diff --git a/doc/index.md b/doc/index.md index e799782e41..04b084811b 100644 --- a/doc/index.md +++ b/doc/index.md @@ -120,7 +120,7 @@ Prepare data for machine learning algorithms. Transform data from one space to another. - * [`PCA`](pca.md): principal components analysis + * [`PCA`](user/methods/pca.md): principal components analysis ### Modeling utilities diff --git a/doc/user/methods/pca.md b/doc/user/methods/pca.md index b01a8404f8..ccf65048d2 100644 --- a/doc/user/methods/pca.md +++ b/doc/user/methods/pca.md @@ -30,7 +30,7 @@ std::cout << "The transformed data matrix has size " << dataset.n_rows /* 5 */ #### Quick links: * [Constructors](#constructors): create `PCA` objects. - * [`Apply()`](#applying-transformation): apply PCA transformation to data. + * [`Apply()`](#applying-transformations): apply PCA transformation to data. * [Examples](#simple-examples) of simple usage and links to detailed example projects. * [Template parameters](#advanced-functionality-different-decomposition-strategies) @@ -38,9 +38,10 @@ std::cout << "The transformed data matrix has size " << dataset.n_rows /* 5 */ #### See also: - * [`RADICAL`](radical.md): independent components analysis - * [mlpack preprocessing utilities](index.md#preprocessing-utilities) - * [mlpack transformations](index.md#transformations) + + + * [mlpack preprocessing utilities](../../index.md#preprocessing-utilities) + * [mlpack transformations](../../index.md#transformations) * [Principal component analysis on Wikipedia](https://en.wikipedia.org/wiki/Principal_component_analysis) ### Constructors From e8792de3bc77fd90d34b6e9afa23e898d6770fe3 Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Wed, 3 Apr 2024 13:57:55 -0400 Subject: [PATCH 08/11] Don't fail on 301s (bans for crawling). --- scripts/build-docs.sh | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/scripts/build-docs.sh b/scripts/build-docs.sh index 49b6360164..104e59310c 100755 --- a/scripts/build-docs.sh +++ b/scripts/build-docs.sh @@ -320,7 +320,7 @@ do echo "Checking links in $f..."; # To run checklink we have to strip out some perl stderr warnings... - checklink -qs --follow-file-links --suppress-broken 405 "$f" 2>&1 | + checklink -qs --follow-file-links --suppress-broken 405 --suppress-broken 301 "$f" 2>&1 | grep -v 'Use of uninitialized value' > checklink_out; if [ -s checklink_out ]; then From fb6baba77982f6f3d5196482f4705eba17451549 Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Wed, 3 Apr 2024 13:58:02 -0400 Subject: [PATCH 09/11] Update Weka links. --- doc/user/load_save.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/doc/user/load_save.md b/doc/user/load_save.md index 9785c968e9..c7fa95627c 100644 --- a/doc/user/load_save.md +++ b/doc/user/load_save.md @@ -157,7 +157,7 @@ With a `data::DatasetInfo` object, categorical data can be loaded: contents of the file: - `.csv`, `.tsv`, or `.txt` for CSV/TSV (tab-separated)/ASCII (space-separated) - - `.arff` for [ARFF](https://www.cs.waikato.ac.nz/~ml/weka/arff.html) + - `.arff` for [ARFF](https://ml.cms.waikato.ac.nz/weka/arff.html) * `matrix` is an `arma::mat&`, `arma::Mat&`, or similar (e.g., a reference to an Armadillo object that data will be loaded into or saved @@ -723,7 +723,7 @@ The format of mixed categorical data is detected automatically based on the file extension and inspecting the file contents: - `.csv`, `.txt`, or `.tsv` indicates CSV/TSV/ASCII format - - `.arff` indicates [ARFF](https://www.cs.waikato.ac.nz/~ml/weka/arff.html) + - `.arff` indicates [ARFF](https://ml.cms.waikato.ac.nz/weka/arff.html) --- From 58d5e0250bf3c2fa97c5dd92fc857d28df76bc2b Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Wed, 3 Apr 2024 14:00:02 -0400 Subject: [PATCH 10/11] Update HISTORY.md. --- HISTORY.md | 2 ++ 1 file changed, 2 insertions(+) diff --git a/HISTORY.md b/HISTORY.md index ce4ec825f0..d3fe453aa2 100644 --- a/HISTORY.md +++ b/HISTORY.md @@ -21,6 +21,8 @@ * Fix floating-point accuracy issue for decision trees that sometimes caused crashes (#3595). + * Allow PCA to take different matrix types (#3677). + ### mlpack 4.3.0 ###### 2023-11-27 * Fix include ordering issue for `LinearRegression` (#3541). From 787e14b41903a631869b8d51990028b2db2a6795 Mon Sep 17 00:00:00 2001 From: Ryan Curtin Date: Mon, 8 Apr 2024 10:40:15 -0400 Subject: [PATCH 11/11] Fix indentation. --- .../block_krylov_svd/randomized_block_krylov_svd_impl.hpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp b/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp index 7a8000ebc0..567dd2ccbf 100644 --- a/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp +++ b/src/mlpack/methods/block_krylov_svd/randomized_block_krylov_svd_impl.hpp @@ -60,8 +60,8 @@ inline void RandomizedBlockKrylovSVD::Apply(const InMatType& data, { // The block size cannot be greater than the number of points in the // dataset or the dimensionality of the dataset. - blockSize = std::min((size_t) data.n_rows, - std::min((size_t) data.n_cols, rank + 10)); + blockSize = std::min((size_t) data.n_rows, std::min((size_t) data.n_cols, + rank + 10)); } // Random block initialization.