From 3bf42f7e6d10bc2999b5b9e8acff79ed9917c42c Mon Sep 17 00:00:00 2001 From: Mauricio 'Pacha' Vargas Sepulveda Date: Wed, 15 Jan 2025 10:03:16 +0000 Subject: [PATCH] 0.4.3 in progress --- NEWS.md | 5 +- cpp11armadillotest/R/cpp11.R | 32 ++ .../src/08_official_documentation_adapted.cpp | 136 ++++++- cpp11armadillotest/src/cpp11.cpp | 64 +++ .../test-official-documentation-adapted.R | 106 +++++ inst/include/wrappers/matrices.hpp | 5 + inst/include/wrappers/vectors.hpp | 5 + .../functions-of-vector-matrices-cubes.Rmd | 10 +- vignettes/statistics-and-clustering.Rmd | 380 ++++++++++++++++++ 9 files changed, 735 insertions(+), 8 deletions(-) diff --git a/NEWS.md b/NEWS.md index 275474d..3a0d334 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,8 +1,11 @@ # cpp11armadillo 0.4.3 * Expanded documentation. -* Added a new vignette "Functions of vectors, matrices, and cubes". +* Added a new vignettes: + * Functions of vectors, matrices, and cubes. + * Statistics and clustering. * New package template. +* Provides `as_mat()` and `as_col()` as wrappers for `as_Mat()` and `as_Col()`. # cpp11armadillo 0.4.2 diff --git a/cpp11armadillotest/R/cpp11.R b/cpp11armadillotest/R/cpp11.R index 6b4067a..b584220 100644 --- a/cpp11armadillotest/R/cpp11.R +++ b/cpp11armadillotest/R/cpp11.R @@ -780,6 +780,38 @@ trig1_ <- function(n) { .Call(`_cpp11armadillotest_trig1_`, n) } +mean1_ <- function(X, Y) { + .Call(`_cpp11armadillotest_mean1_`, X, Y) +} + +median1_ <- function(X, Y) { + .Call(`_cpp11armadillotest_median1_`, X, Y) +} + +stddev1_ <- function(X, Y) { + .Call(`_cpp11armadillotest_stddev1_`, X, Y) +} + +var1_ <- function(X, Y) { + .Call(`_cpp11armadillotest_var1_`, X, Y) +} + +range1_ <- function(X, Y) { + .Call(`_cpp11armadillotest_range1_`, X, Y) +} + +cov1_ <- function(X, Y) { + .Call(`_cpp11armadillotest_cov1_`, X, Y) +} + +cor1_ <- function(X, Y) { + .Call(`_cpp11armadillotest_cor1_`, X, Y) +} + +hist1_ <- function(n) { + .Call(`_cpp11armadillotest_hist1_`, n) +} + ols_ <- function(x, y) { .Call(`_cpp11armadillotest_ols_`, x, y) } diff --git a/cpp11armadillotest/src/08_official_documentation_adapted.cpp b/cpp11armadillotest/src/08_official_documentation_adapted.cpp index a7753f3..f35e2c1 100644 --- a/cpp11armadillotest/src/08_official_documentation_adapted.cpp +++ b/cpp11armadillotest/src/08_official_documentation_adapted.cpp @@ -2119,13 +2119,16 @@ return res; } -// TODO: check later why vec + cx_vec reurns an error [[cpp11::register]] list roots1_(const int& n) { // y = p_1*x^n + p_2*x^(n-1) + ... + p_(n-1)*x + p_n // p_1, ..., p_n are random numbers - mat A(n, 1, fill::randu); - cx_mat B = roots(A); - list res = as_complex_matrix(B); + vec y(n, 1, fill::randu); + + // note that mat and cx_mat operate directly + // but vec and cx_vec require conv_to<...>::from() + cx_vec z = roots(conv_to::from(y)); + + list res = as_complex_doubles(z); return res; } @@ -2386,3 +2389,128 @@ return res; } + +[[cpp11::register]] list mean1_(const doubles_matrix<>& X, const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + // create a cube with 3 copies of B + random noise + cube C(B.n_rows, B.n_cols, 3); + C.slice(0) = B + 0.1 * randn(B.n_rows, B.n_cols); + C.slice(1) = B + 0.2 * randn(B.n_rows, B.n_cols); + C.slice(2) = B + 0.3 * randn(B.n_rows, B.n_cols); + + vec D = mean(A).t(); + vec E = mean(A, 1); + vec F = mean(mean(B, 1), 1); + + writable::list res(3); + res[0] = as_doubles(D); + res[1] = as_doubles(E); + res[2] = as_doubles(F); + + return res; +} + +[[cpp11::register]] list median1_(const doubles_matrix<>& X, const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + vec C = median(A).t(); + vec D = median(A, 1); + vec E = median(median(B, 1), 1); + + writable::list res(3); + res[0] = as_doubles(C); + res[1] = as_doubles(D); + res[2] = as_doubles(E); + + return res; +} + +[[cpp11::register]] list stddev1_(const doubles_matrix<>& X, const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + vec C = stddev(A).t(); + vec D = stddev(A, 1).t(); + vec E = stddev(A, 1, 1); + + writable::list res(3); + res[0] = as_doubles(C); + res[1] = as_doubles(D); + res[2] = as_doubles(E); + + return res; +} + +[[cpp11::register]] list var1_(const doubles_matrix<>& X, const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + vec C = var(A).t(); + vec D = var(A, 1).t(); + vec E = var(A, 1, 1); + + writable::list res(3); + res[0] = as_doubles(C); + res[1] = as_doubles(D); + res[2] = as_doubles(E); + + return res; +} + +[[cpp11::register]] list range1_(const doubles_matrix<>& X, const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + vec C = range(A).t(); + vec D = range(A, 1); + + writable::list res(2); + res[0] = as_doubles(C); + res[1] = as_doubles(D); + + return res; +} + +[[cpp11::register]] list cov1_(const doubles_matrix<>& X, const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + mat C = cov(A, B); + mat D = cov(A, B, 1); + + writable::list res(2); + res[0] = as_doubles_matrix(C); + res[1] = as_doubles_matrix(D); + + return res; +} + +[[cpp11::register]] list cor1_(const doubles_matrix<>& X, const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + mat C = cor(A, B); + mat D = cor(A, B, 1); + + writable::list res(2); + res[0] = as_doubles_matrix(C); + res[1] = as_doubles_matrix(D); + + return res; +} + +[[cpp11::register]] list hist1_(const int& n) { + vec A = randu(n); + + uvec h1 = hist(A, 11); + uvec h2 = hist(A, linspace(-2, 2, 11)); + + writable::list res(2); + res[0] = as_integers(h1); + res[1] = as_integers(h2); + + return res; +} diff --git a/cpp11armadillotest/src/cpp11.cpp b/cpp11armadillotest/src/cpp11.cpp index cc611d9..c412be2 100644 --- a/cpp11armadillotest/src/cpp11.cpp +++ b/cpp11armadillotest/src/cpp11.cpp @@ -1370,6 +1370,62 @@ extern "C" SEXP _cpp11armadillotest_trig1_(SEXP n) { return cpp11::as_sexp(trig1_(cpp11::as_cpp>(n))); END_CPP11 } +// 08_official_documentation_adapted.cpp +list mean1_(const doubles_matrix<>& X, const doubles_matrix<>& Y); +extern "C" SEXP _cpp11armadillotest_mean1_(SEXP X, SEXP Y) { + BEGIN_CPP11 + return cpp11::as_sexp(mean1_(cpp11::as_cpp&>>(X), cpp11::as_cpp&>>(Y))); + END_CPP11 +} +// 08_official_documentation_adapted.cpp +list median1_(const doubles_matrix<>& X, const doubles_matrix<>& Y); +extern "C" SEXP _cpp11armadillotest_median1_(SEXP X, SEXP Y) { + BEGIN_CPP11 + return cpp11::as_sexp(median1_(cpp11::as_cpp&>>(X), cpp11::as_cpp&>>(Y))); + END_CPP11 +} +// 08_official_documentation_adapted.cpp +list stddev1_(const doubles_matrix<>& X, const doubles_matrix<>& Y); +extern "C" SEXP _cpp11armadillotest_stddev1_(SEXP X, SEXP Y) { + BEGIN_CPP11 + return cpp11::as_sexp(stddev1_(cpp11::as_cpp&>>(X), cpp11::as_cpp&>>(Y))); + END_CPP11 +} +// 08_official_documentation_adapted.cpp +list var1_(const doubles_matrix<>& X, const doubles_matrix<>& Y); +extern "C" SEXP _cpp11armadillotest_var1_(SEXP X, SEXP Y) { + BEGIN_CPP11 + return cpp11::as_sexp(var1_(cpp11::as_cpp&>>(X), cpp11::as_cpp&>>(Y))); + END_CPP11 +} +// 08_official_documentation_adapted.cpp +list range1_(const doubles_matrix<>& X, const doubles_matrix<>& Y); +extern "C" SEXP _cpp11armadillotest_range1_(SEXP X, SEXP Y) { + BEGIN_CPP11 + return cpp11::as_sexp(range1_(cpp11::as_cpp&>>(X), cpp11::as_cpp&>>(Y))); + END_CPP11 +} +// 08_official_documentation_adapted.cpp +list cov1_(const doubles_matrix<>& X, const doubles_matrix<>& Y); +extern "C" SEXP _cpp11armadillotest_cov1_(SEXP X, SEXP Y) { + BEGIN_CPP11 + return cpp11::as_sexp(cov1_(cpp11::as_cpp&>>(X), cpp11::as_cpp&>>(Y))); + END_CPP11 +} +// 08_official_documentation_adapted.cpp +list cor1_(const doubles_matrix<>& X, const doubles_matrix<>& Y); +extern "C" SEXP _cpp11armadillotest_cor1_(SEXP X, SEXP Y) { + BEGIN_CPP11 + return cpp11::as_sexp(cor1_(cpp11::as_cpp&>>(X), cpp11::as_cpp&>>(Y))); + END_CPP11 +} +// 08_official_documentation_adapted.cpp +list hist1_(const int& n); +extern "C" SEXP _cpp11armadillotest_hist1_(SEXP n) { + BEGIN_CPP11 + return cpp11::as_sexp(hist1_(cpp11::as_cpp>(n))); + END_CPP11 +} // 09_regression.cpp doubles ols_(const doubles_matrix<>& x, const doubles& y); extern "C" SEXP _cpp11armadillotest_ols_(SEXP x, SEXP y) { @@ -1414,6 +1470,8 @@ static const R_CallMethodDef CallEntries[] = { {"_cpp11armadillotest_conj1_", (DL_FUNC) &_cpp11armadillotest_conj1_, 1}, {"_cpp11armadillotest_conv_to1_", (DL_FUNC) &_cpp11armadillotest_conv_to1_, 1}, {"_cpp11armadillotest_copy_size1_", (DL_FUNC) &_cpp11armadillotest_copy_size1_, 1}, + {"_cpp11armadillotest_cor1_", (DL_FUNC) &_cpp11armadillotest_cor1_, 2}, + {"_cpp11armadillotest_cov1_", (DL_FUNC) &_cpp11armadillotest_cov1_, 2}, {"_cpp11armadillotest_cross1_", (DL_FUNC) &_cpp11armadillotest_cross1_, 1}, {"_cpp11armadillotest_cube1_", (DL_FUNC) &_cpp11armadillotest_cube1_, 2}, {"_cpp11armadillotest_cumprod1_", (DL_FUNC) &_cpp11armadillotest_cumprod1_, 1}, @@ -1450,6 +1508,7 @@ static const R_CallMethodDef CallEntries[] = { {"_cpp11armadillotest_for_each1_", (DL_FUNC) &_cpp11armadillotest_for_each1_, 1}, {"_cpp11armadillotest_has_inf1_", (DL_FUNC) &_cpp11armadillotest_has_inf1_, 1}, {"_cpp11armadillotest_has_nan1_", (DL_FUNC) &_cpp11armadillotest_has_nan1_, 1}, + {"_cpp11armadillotest_hist1_", (DL_FUNC) &_cpp11armadillotest_hist1_, 1}, {"_cpp11armadillotest_imag1_", (DL_FUNC) &_cpp11armadillotest_imag1_, 1}, {"_cpp11armadillotest_imbue1_", (DL_FUNC) &_cpp11armadillotest_imbue1_, 1}, {"_cpp11armadillotest_imbue2_", (DL_FUNC) &_cpp11armadillotest_imbue2_, 1}, @@ -1493,6 +1552,8 @@ static const R_CallMethodDef CallEntries[] = { {"_cpp11armadillotest_matrix2_", (DL_FUNC) &_cpp11armadillotest_matrix2_, 1}, {"_cpp11armadillotest_max1_", (DL_FUNC) &_cpp11armadillotest_max1_, 1}, {"_cpp11armadillotest_maxmin1_", (DL_FUNC) &_cpp11armadillotest_maxmin1_, 1}, + {"_cpp11armadillotest_mean1_", (DL_FUNC) &_cpp11armadillotest_mean1_, 2}, + {"_cpp11armadillotest_median1_", (DL_FUNC) &_cpp11armadillotest_median1_, 2}, {"_cpp11armadillotest_memptr1_", (DL_FUNC) &_cpp11armadillotest_memptr1_, 1}, {"_cpp11armadillotest_misc1_", (DL_FUNC) &_cpp11armadillotest_misc1_, 1}, {"_cpp11armadillotest_nonzeros1_", (DL_FUNC) &_cpp11armadillotest_nonzeros1_, 1}, @@ -1519,6 +1580,7 @@ static const R_CallMethodDef CallEntries[] = { {"_cpp11armadillotest_randu1_", (DL_FUNC) &_cpp11armadillotest_randu1_, 1}, {"_cpp11armadillotest_randu2_", (DL_FUNC) &_cpp11armadillotest_randu2_, 1}, {"_cpp11armadillotest_randu3_", (DL_FUNC) &_cpp11armadillotest_randu3_, 1}, + {"_cpp11armadillotest_range1_", (DL_FUNC) &_cpp11armadillotest_range1_, 2}, {"_cpp11armadillotest_rank1_", (DL_FUNC) &_cpp11armadillotest_rank1_, 1}, {"_cpp11armadillotest_rcond1_", (DL_FUNC) &_cpp11armadillotest_rcond1_, 1}, {"_cpp11armadillotest_regspace1_", (DL_FUNC) &_cpp11armadillotest_regspace1_, 1}, @@ -1553,6 +1615,7 @@ static const R_CallMethodDef CallEntries[] = { {"_cpp11armadillotest_sprandu1_", (DL_FUNC) &_cpp11armadillotest_sprandu1_, 1}, {"_cpp11armadillotest_sqrtmat1_", (DL_FUNC) &_cpp11armadillotest_sqrtmat1_, 1}, {"_cpp11armadillotest_sqrtmat_sympd1_", (DL_FUNC) &_cpp11armadillotest_sqrtmat_sympd1_, 1}, + {"_cpp11armadillotest_stddev1_", (DL_FUNC) &_cpp11armadillotest_stddev1_, 2}, {"_cpp11armadillotest_sub2ind1_", (DL_FUNC) &_cpp11armadillotest_sub2ind1_, 1}, {"_cpp11armadillotest_subview1_", (DL_FUNC) &_cpp11armadillotest_subview1_, 1}, {"_cpp11armadillotest_subview2_", (DL_FUNC) &_cpp11armadillotest_subview2_, 1}, @@ -1580,6 +1643,7 @@ static const R_CallMethodDef CallEntries[] = { {"_cpp11armadillotest_typedef_SpMat_int", (DL_FUNC) &_cpp11armadillotest_typedef_SpMat_int, 1}, {"_cpp11armadillotest_typedef_uvec", (DL_FUNC) &_cpp11armadillotest_typedef_uvec, 1}, {"_cpp11armadillotest_unique1_", (DL_FUNC) &_cpp11armadillotest_unique1_, 1}, + {"_cpp11armadillotest_var1_", (DL_FUNC) &_cpp11armadillotest_var1_, 2}, {"_cpp11armadillotest_vecnorm1_", (DL_FUNC) &_cpp11armadillotest_vecnorm1_, 1}, {"_cpp11armadillotest_vectorise1_", (DL_FUNC) &_cpp11armadillotest_vectorise1_, 1}, {"_cpp11armadillotest_zeros1_", (DL_FUNC) &_cpp11armadillotest_zeros1_, 1}, diff --git a/cpp11armadillotest/tests/testthat/test-official-documentation-adapted.R b/cpp11armadillotest/tests/testthat/test-official-documentation-adapted.R index 0e7caf3..39d0a12 100644 --- a/cpp11armadillotest/tests/testthat/test-official-documentation-adapted.R +++ b/cpp11armadillotest/tests/testthat/test-official-documentation-adapted.R @@ -558,76 +558,182 @@ test_that("examples derived from official documentation", { res128 <- inplace_strans1_(2) expect_type(res128, "list") + expect_equal(length(res128), 2) res129 <- intersect1_(5) expect_type(res129, "integer") + expect_equal(res129, 2:5) res130 <- join_rows1_(2) expect_type(res130, "list") + expect_equal(length(res130), 5) res131 <- join_cubes1_(2) expect_type(res131, "list") + expect_equal(length(res131), 7) res132 <- kron1_(2) expect_type(res132, "double") + expect_equal(dim(res132), c(6, 6)) res133 <- log_det1_(2) expect_type(res133, "list") + expect_equal(length(res133), 2) res134 <- log_det_sympd1_(2) expect_type(res134, "list") + expect_equal(length(res134), 2) res135 <- logmat1_(2) expect_type(res135, "list") + expect_equal(length(res135), 2) res136 <- logmat_sympd1_(2) expect_type(res136, "list") + expect_equal(length(res136), 2) res137 <- max1_(2) expect_type(res137, "list") + expect_equal(length(res137), 4) res138 <- nonzeros1_(2) expect_type(res138, "double") res139 <- norm1_(2) expect_type(res139, "double") + expect_equal(length(res139), 5) res140 <- norm2est1_(2) expect_type(res140, "double") + expect_equal(length(res140), 1) res141 <- normalise1_(2) expect_type(res141, "list") + expect_equal(length(res141), 2) res142 <- pow1_(2) expect_type(res142, "list") + expect_equal(length(res142), 2) res143 <- powmat1_(2) expect_type(res143, "list") + expect_equal(length(res143), 2) res144 <- prod1_(2) expect_type(res144, "list") + expect_equal(length(res144), 2) res145 <- rank1_(2) expect_type(res145, "list") + expect_equal(length(res145), 2) res146 <- rcond1_(2) expect_type(res146, "double") + expect_equal(length(res146), 1) res147 <- repelem1_(2) expect_type(res147, "list") + expect_equal(length(res147), 2) res148 <- repmat1_(2) expect_type(res148, "list") + expect_equal(length(res148), 2) res149 <- reshape2_(2) expect_type(res149, "list") + expect_equal(length(res149), 3) res150 <- resize2_(2) expect_type(res150, "list") + expect_equal(length(res150), 3) res151 <- reverse1_(2) expect_type(res151, "list") + expect_equal(length(res151), 3) res152 <- roots1_(2) expect_type(res152, "list") + expect_equal(length(res152), 2) + + res153 <- shift1_(2) + expect_type(res153, "list") + expect_equal(length(res153), 3) + + res154 <- shuffle1_(2) + expect_type(res154, "list") + expect_equal(length(res154), 2) + + res155 <- size1_(2) + expect_type(res155, "list") + expect_equal(length(res155), 7) + + res156 <- sort1_(2) + expect_type(res156, "list") + expect_equal(length(res156), 5) + + res157 <- sort_index1_(2) + expect_type(res157, "list") + expect_equal(length(res157), 3) + + res158 <- sqrtmat1_(2) + expect_type(res158, "list") + expect_equal(length(res158), 4) + + res159 <- sqrtmat_sympd1_(2) + expect_type(res159, "double") + expect_equal(dim(res159), c(2, 2)) + + res160 <- sum2_(2) + expect_type(res160, "list") + expect_equal(length(res160), 3) + + res161 <- sub2ind1_(2) + expect_type(res161, "integer") + expect_equal(length(res161), 1) + + res162 <- trace1_(2) + expect_type(res162, "double") + expect_equal(length(res162), 1) + + res163 <- trans1_(2) + expect_type(res163, "list") + expect_equal(length(res163), 2) + + res164 <- trapz1_(2) + expect_type(res164, "double") + expect_equal(length(res164), 1) + + set.seed(123) + X <- matrix(rnorm(4), nrow = 2, ncol = 2) + Y <- matrix(rnorm(4), nrow = 2, ncol = 2) + + res165 <- mean1_(X, Y) + res166 <- median1_(X, Y) + res167 <- stddev1_(X, Y) + res168 <- var1_(X, Y) + res169 <- range1_(X, Y) + res170 <- cov1_(X,Y) + res171 <- cor1_(X,Y) + + expect_type(res165, "list") + expect_equal(length(res165), 3) + + expect_type(res166, "list") + expect_equal(length(res166), 3) + expect_equal(res165, res166) + + expect_type(res167, "list") + expect_equal(length(res167), 3) + + expect_type(res168, "list") + expect_equal(length(res168), 3) + + expect_type(res169, "list") + expect_equal(length(res169), 2) + + expect_type(res170, "list") + expect_equal(length(res170), 2) + + expect_type(res171, "list") + expect_equal(length(res171), 2) }) diff --git a/inst/include/wrappers/matrices.hpp b/inst/include/wrappers/matrices.hpp index c45d464..59a52e4 100644 --- a/inst/include/wrappers/matrices.hpp +++ b/inst/include/wrappers/matrices.hpp @@ -16,6 +16,11 @@ inline Mat as_Mat(const T& x) { throw std::runtime_error("Cannot convert to Mat"); } +// armadillo 0.4.3 +// as_mat = as_Mat +template +inline Mat as_Mat(const Mat& x) { return x; } + template inline Mat dblint_matrix_to_Mat_(const U& x) { const int n = x.nrow(); diff --git a/inst/include/wrappers/vectors.hpp b/inst/include/wrappers/vectors.hpp index df1a442..d2eff98 100644 --- a/inst/include/wrappers/vectors.hpp +++ b/inst/include/wrappers/vectors.hpp @@ -20,6 +20,11 @@ inline Col as_Col(const T& x) { throw std::runtime_error("Cannot convert to Col"); } +// armadillo 0.4.3 +// as_col = as_Col +template +inline Col as_Col(const Col& x) { return x; } + template inline Col as_Col_(const U& x) { const size_t n = x.size(); diff --git a/vignettes/functions-of-vector-matrices-cubes.Rmd b/vignettes/functions-of-vector-matrices-cubes.Rmd index f49b255..98b33d3 100644 --- a/vignettes/functions-of-vector-matrices-cubes.Rmd +++ b/vignettes/functions-of-vector-matrices-cubes.Rmd @@ -1861,9 +1861,13 @@ roots(Y, X) // store the roots in Y and return true if successful [[cpp11::register]] list roots1_(const int& n) { // y = p_1*x^n + p_2*x^(n-1) + ... + p_(n-1)*x + p_n // p_1, ..., p_n are random numbers - mat A(n, 1, fill::randu); - cx_mat B = roots(A); - list res = as_complex_matrix(B); + vec y(n, 1, fill::randu); + + // note that mat and cx_mat operate directly + // but vec and cx_vec require conv_to<...>::from() + cx_vec z = roots(conv_to::from(y)); + + list res = as_complex_doubles(z); return res; } ``` diff --git a/vignettes/statistics-and-clustering.Rmd b/vignettes/statistics-and-clustering.Rmd index 51dd909..266f23f 100644 --- a/vignettes/statistics-and-clustering.Rmd +++ b/vignettes/statistics-and-clustering.Rmd @@ -17,3 +17,383 @@ knitr::opts_chunk$set( This vignette is adapted from the official Armadillo [documentation](https://arma.sourceforge.net/docs.html). + +| Function | Description | +|--------------------|---------------------------------------------------------| +| `mean` | Mean {#mean} | +| `median` | Median {#median} | +| `stddev` | Standard deviation {#stddev} | +| `var` | Variance {#var} | +| `range` | Range {#range} | +| `cov` | Covariance {#cov} | +| `cor` | Correlation {#cor} | +| `hist` | Histogram of counts {#hist} | +| `histc` | Histogram of counts with user specified edges {#histc} | +| `quantile` | Quantiles of a dataset {#quantile} | +| `princomp` | Principal component analysis (PCA) {#princomp} | +| `normpdf` | Probability density function of normal distribution {#normpdf} | +| `log_normpdf` | Logarithm version of probability density function of normal distribution {#log_normpdf} | +| `normcdf` | Cumulative distribution function of normal distribution {#normcdf} | +| `mvnrnd` | Random vectors from multivariate normal distribution {#mvnrnd} | +| `chi2rnd` | Random numbers from chi-squared distribution {#chi2rnd} | +| `wishrnd` | Random matrix from Wishart distribution {#wishrnd} | +| `iwishrnd` | Random matrix from inverse Wishart distribution {#iwishrnd} | +| `running_stat` | Running statistics of scalars (one dimensional process/signal) {#running_stat} | +| `running_stat_vec` | Running statistics of vectors (multi-dimensional process/signal) {#running_stat_vec} | +| `kmeans` | Cluster data into disjoint sets {#kmeans} | +| `gmm_diag` / `gmm_full` | Probabilistic clustering and likelihood calculation via mixture of Gaussians {#gmm_diag-gmm_full} | + +# Mean {#mean} + +The `mean` function computes the mean of a vector, matrix, or cube. For a vector +argument, the mean is calculated using all the elements of the vector. For a +matrix argument, the mean is calculated for each column by default (`dim = 0`), +or for each row (`dim = 1`). For a cube argument, the mean is calculated along +the specified dimension (same as for matrices plus `dim = 2` for slices). + +Usage: + +```cpp +mean(V) + +mean(M) +mean(M, dim) + +mean(Q) +mean(Q, dim) +``` + +## Examples + +```cpp +[[cpp11::register]] list mean1_(const doubles_matrix<>& X, + const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + // create a cube with 3 copies of B + random noise + cube C(B.n_rows, B.n_cols, 3); + C.slice(0) = B + 0.1 * randn(B.n_rows, B.n_cols); + C.slice(1) = B + 0.2 * randn(B.n_rows, B.n_cols); + C.slice(2) = B + 0.3 * randn(B.n_rows, B.n_cols); + + vec D = mean(A).t(); + vec E = mean(A, 1); + vec F = mean(mean(B, 1), 1); + + writable::list res(3); + res[0] = as_doubles(D); + res[1] = as_doubles(E); + res[2] = as_doubles(F); + + return res; +} +``` + +# Median {#median} + +The `median` function computes the median of a vector or matrix. For a vector +argument, the median is calculated using all the elements of the vector. For a +matrix argument, the median is calculated for each column by default +(`dim = 0`), or for each row (`dim = 1`). + +Usage: + +```cpp +median(V) + +median(M) +median(M, dim) +``` + +## Examples + +```cpp +[[cpp11::register]] list median1_(const doubles_matrix<>& X, + const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + vec C = median(A).t(); + vec D = median(A, 1); + vec E = median(median(B, 1), 1); + + writable::list res(3); + res[0] = as_doubles(C); + res[1] = as_doubles(D); + res[2] = as_doubles(E); + + return res; +} +``` + +# Standard deviation {#stddev} + +The `stddev` function computes the standard deviation of a vector or matrix. +For a vector argument, the standard deviation is calculated using all the +elements of the vector. For a matrix argument, the standard deviation is +calculated for each column by default (`dim = 0`), or for each row (`dim = 1`). + +The `norm_type` argument is optional; by default `norm_type = 0` is used. The +`norm_type` argument controls the type of normalization used, with `N` denoting +the number of observations: + +- for `norm_type = 0`, normalization is done using `N-1`, providing the best + unbiased estimation of the standard deviation (if the observations are from a + normal distribution). +- for `norm_type = 1`, normalization is done using `N`, which provides the + second moment of the observations about their mean. + +Usage: + +```cpp +stddev(V) +stddev(V, norm_type) + +stddev(M) +stddev(M, norm_type) +stddev(M, norm_type, dim) +``` + +## Examples + +```cpp +[[cpp11::register]] list stddev1_(const doubles_matrix<>& X, + const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + vec C = stddev(A).t(); + vec D = stddev(A, 1).t(); + vec E = stddev(A, 1, 1); + + writable::list res(3); + res[0] = as_doubles(C); + res[1] = as_doubles(D); + res[2] = as_doubles(E); + + return res; +} +``` + +# Variance {#var} + +The `var` function computes the variance of a vector or matrix. +For a vector argument, the variance is calculated using all the +elements of the vector. For a matrix argument, the variance is +calculated for each column by default (`dim = 0`), or for each row (`dim = 1`). + +The `norm_type` argument is optional; by default `norm_type = 0` is used. The +`norm_type` argument controls the type of normalization used, with `N` denoting +the number of observations: + +- for `norm_type = 0`, normalization is done using `N-1`, providing the best + unbiased estimation of the standard deviation (if the observations are from a + normal distribution). +- for `norm_type = 1`, normalization is done using `N`, which provides the + second moment of the observations about their mean. + +Usage: + +```cpp +var(V) +var(V, norm_type) + +var(M) +var(M, norm_type) +var(M, norm_type, dim) +``` + +## Examples + +```cpp +[[cpp11::register]] list var1_(const doubles_matrix<>& X, + const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + vec C = var(A).t(); + vec D = var(A, 1).t(); + vec E = var(A, 1, 1); + + writable::list res(3); + res[0] = as_doubles(C); + res[1] = as_doubles(D); + res[2] = as_doubles(E); + + return res; +} +``` + +# Range {#range} + +The `range` function computes the range of a vector or matrix. For a vector +argument, the range is calculated using all the elements of the vector. For a +matrix argument, the range is calculated for each column by default (`dim = 0`), +or for each row (`dim = 1`). + +Usage: + +```cpp +range(V) + +range(M) +range(M, dim) +``` + +## Examples + +```cpp +[[cpp11::register]] list range1_(const doubles_matrix<>& X, + const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + vec C = range(A).t(); + vec D = range(A, 1); + + writable::list res(2); + res[0] = as_doubles(C); + res[1] = as_doubles(D); + + return res; +} +``` + +# Covariance {#cov} + +The `cov` function computes the covariance between two matrices or vectors. If +each row of `X` and `Y` is an observation and each column is a variable, the +`(i,j)`-th entry of `cov(X,Y)` is the covariance between the `i`-th variable in +`X` and the `j`-th variable in `Y`. + +For two matrix arguments `X` and `Y`, the `cov(X,Y)` function computes the +covariance between the two matrices. For vector arguments, the type of vector is +ignored and each element in the vector is treated as an observation. The +`cov(X)` function is equivalent to `cov(X, X)`. + +The `norm_type` argument is optional; by default `norm_type = 0` is used. The +`norm_type` argument controls the type of normalization used, with `N` denoting +the number of observations: + +- for `norm_type = 0`, normalization is done using `N-1`, providing the best + unbiased estimation of the covariance matrix (if the observations are from a + normal distribution). +- for `norm_type = 1`, normalization is done using `N`, which provides the + second moment matrix of the observations about their mean. + +Usage: + +```cpp +cov(X, Y, norm_type) + +cov(X) +cov(X, norm_type) +``` + +## Examples + +```cpp +[[cpp11::register]] list cov1_(const doubles_matrix<>& X, + const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + mat C = cov(A, B); + mat D = cov(A, B, 1); + + writable::list res(2); + res[0] = as_doubles_matrix(C); + res[1] = as_doubles_matrix(D); + + return res; +} +``` + +# Correlation {#cor} + +The `cor` function computes the correlation coefficient between two matrices or +vectors. If each row of `X` and `Y` is an observation and each column is a +variable, the `(i,j)`-th entry of `cor(X,Y)` is the correlation coefficient +between the `i`-th variable in `X` and the `j`-th variable in `Y`. + +For two matrix arguments `X` and `Y`, the `cor(X,Y)` function computes the +correlation coefficient between the two matrices. For vector arguments, the type +of vector is ignored and each element in the vector is treated as an observation. + +The `norm_type` argument is optional; by default `norm_type = 0` is used. The +`norm_type` argument controls the type of normalization used, with `N` denoting +the number of observations: + +- for `norm_type = 0`, normalization is done using `N-1`. +- for `norm_type = 1`, normalization is done using `N`. + +Usage: + +```cpp +cor(X, Y) +cor(X, Y, norm_type) + +cor(X) +cor(X, norm_type) +``` + +## Examples + +```cpp +[[cpp11::register]] list cor1_(const doubles_matrix<>& X, + const doubles_matrix<>& Y) { + mat A = as_Mat(X); + mat B = as_Mat(Y); + + mat C = cor(A, B); + mat D = cor(A, B, 1); + + writable::list res(2); + res[0] = as_doubles_matrix(C); + res[1] = as_doubles_matrix(D); + + return res; +} +``` + +# Histogram {#hist} + +The `hist` function computes a histogram of counts for a vector or matrix. For a +vector argument, the histogram is calculated using all the elements of the +vector. For a matrix argument, the histogram is calculated for each column by +default (`dim = 0`), or for each row (`dim = 1`). + +The bin centers can be automatically determined from the data, with the number +of bins specified via `n_bins` (the default is 10). The range of the bins is +determined by the range of the data. The bin centers can also be explicitly +specified via the `centers` vector, which must contain monotonically increasing +values. + +Usage: + +```cpp +hist(V) +hist(V, n_bins) +hist(V, centers) + +hist(M, centers) +hist(M, centers, dim) +``` + +## Examples + +```cpp +[[cpp11::register]] list hist1_(const int& n) { + vec A = randu(n); + + uvec h1 = hist(A, 11); + uvec h2 = hist(A, linspace(-2, 2, 11)); + + writable::list res(2); + res[0] = as_integers(h1); + res[1] = as_integers(h2); + + return res; +} +```