From d0aa4c250142307eb274af1616aa55e951a1ff76 Mon Sep 17 00:00:00 2001 From: David Cortes Date: Sun, 29 Nov 2020 22:17:43 +0200 Subject: [PATCH 1/2] re-add mean centering and bias init --- R/RcppExports.R | 4 +-- R/model_WRMF.R | 83 ++++++++++++++++--------------------------- src/RcppExports.cpp | 14 ++++---- src/SolverAlsWrmf.cpp | 56 +++++++++++++++++------------ 4 files changed, 73 insertions(+), 84 deletions(-) diff --git a/R/RcppExports.R b/R/RcppExports.R index 4a0bddf..37caf00 100644 --- a/R/RcppExports.R +++ b/R/RcppExports.R @@ -90,11 +90,11 @@ als_explicit_float <- function(m_csc_r, X_, Y_, lambda, n_threads, solver, cg_st } initialize_biases_double <- function(m_csc_r, m_csr_r, user_bias, item_bias, lambda, non_negative) { - invisible(.Call(`_rsparse_initialize_biases_double`, m_csc_r, m_csr_r, user_bias, item_bias, lambda, non_negative)) + .Call(`_rsparse_initialize_biases_double`, m_csc_r, m_csr_r, user_bias, item_bias, lambda, non_negative) } initialize_biases_float <- function(m_csc_r, m_csr_r, user_bias, item_bias, lambda, non_negative) { - invisible(.Call(`_rsparse_initialize_biases_float`, m_csc_r, m_csr_r, user_bias, item_bias, lambda, non_negative)) + .Call(`_rsparse_initialize_biases_float`, m_csc_r, m_csr_r, user_bias, item_bias, lambda, non_negative) } deep_copy <- function(x) { diff --git a/R/model_WRMF.R b/R/model_WRMF.R index 9c2df54..cc7840c 100644 --- a/R/model_WRMF.R +++ b/R/model_WRMF.R @@ -164,24 +164,6 @@ WRMF = R6::R6Class( stopifnot(all(c_ui@x >= 0)) } - # if (private$feedback == "explicit" && !private$non_negative) { - # self$global_mean = mean(c_ui@x) - # c_ui@x = c_ui@x - self$global_mean - # } - # if (private$with_bias) { - # c_ui@x = deep_copy(c_ui@x) - # c_ui_orig = deep_copy(c_ui@x) - # } - # else { - # c_ui_orig = numeric(0L) - # } - - # if (private$with_bias) { - # c_iu_orig = deep_copy(c_iu@x) - # } else { - # c_iu_orig = numeric(0L) - # } - # init n_user = nrow(c_ui) n_item = ncol(c_ui) @@ -236,28 +218,34 @@ WRMF = R6::R6Class( stopifnot(ncol(private$U) == ncol(c_iu)) stopifnot(ncol(self$components) == ncol(c_ui)) - # if (private$with_bias) { - # logger$debug("initializing biases") - # if (private$precision == "double") { - # user_bias = numeric(n_user) - # item_bias = numeric(n_item) - # initialize_biases_double(c_ui, c_iu, - # user_bias, - # item_bias, - # private$lambda, - # private$non_negative) - # } else { - # user_bias = float(n_user) - # item_bias = float(n_item) - # initialize_biases_float(c_ui, c_iu, - # user_bias, - # item_bias, - # private$lambda, - # private$non_negative) - # } - # self$components[1L, ] = item_bias - # private$U[private$rank, ] = user_bias - # } + if (private$with_bias) { + logger$debug("initializing biases") + if (private$feedback == "explicit") + c_ui@x = deep_copy(c_ui@x) + if (private$precision == "double") { + user_bias = numeric(n_user) + item_bias = numeric(n_item) + self$global_mean = initialize_biases_double(c_ui, c_iu, + user_bias, + item_bias, + private$lambda, + private$non_negative) + } else { + user_bias = float(n_user) + item_bias = float(n_item) + self$global_mean = initialize_biases_float(c_ui, c_iu, + user_bias, + item_bias, + private$lambda, + private$non_negative) + } + self$components[1L, ] = item_bias + private$U[private$rank, ] = user_bias + } else if (private$feedback == "explicit") { + self$global_mean = mean(c_ui@x) + c_ui@x = c_ui@x - self$global_mean + c_iu@x = c_iu@x - self$global_mean + } logger$info("starting factorization with %d threads", getOption("rsparse_omp_threads", 1L)) @@ -317,6 +305,8 @@ WRMF = R6::R6Class( x = as(x, "CsparseMatrix") x = private$preprocess(x) + if (self$global_mean != 0.) + x@x = x@x - self$global_mean if (private$precision == "double") { res = matrix(0, nrow = private$rank, ncol = nrow(x)) @@ -325,19 +315,6 @@ WRMF = R6::R6Class( } loss = private$solver(t(x), self$components, res, FALSE, private$XtX) - # if (private$feedback == "implicit") { - # loss = private$solver(t(x), self$components, res, FALSE, private$XtX) - # } else{ - # x_use = t(x) - # if (!private$non_negative) - # x_use@x = x_use@x - self$global_mean - # if (private$with_bias) { - # x_orig = deep_copy(x_use@x) - # } else { - # x_orig = numeric(0L) - # } - # loss = private$solver(x_use, self$components, res, FALSE) - # } res = t(res) if (private$precision == "double") diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index ad5c01a..78af4c1 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -322,9 +322,10 @@ BEGIN_RCPP END_RCPP } // initialize_biases_double -void initialize_biases_double(const Rcpp::S4& m_csc_r, const Rcpp::S4& m_csr_r, arma::Col& user_bias, arma::Col& item_bias, double lambda, bool non_negative); +double initialize_biases_double(const Rcpp::S4& m_csc_r, const Rcpp::S4& m_csr_r, arma::Col& user_bias, arma::Col& item_bias, double lambda, bool non_negative); RcppExport SEXP _rsparse_initialize_biases_double(SEXP m_csc_rSEXP, SEXP m_csr_rSEXP, SEXP user_biasSEXP, SEXP item_biasSEXP, SEXP lambdaSEXP, SEXP non_negativeSEXP) { BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; Rcpp::RNGScope rcpp_rngScope_gen; Rcpp::traits::input_parameter< const Rcpp::S4& >::type m_csc_r(m_csc_rSEXP); Rcpp::traits::input_parameter< const Rcpp::S4& >::type m_csr_r(m_csr_rSEXP); @@ -332,14 +333,15 @@ BEGIN_RCPP Rcpp::traits::input_parameter< arma::Col& >::type item_bias(item_biasSEXP); Rcpp::traits::input_parameter< double >::type lambda(lambdaSEXP); Rcpp::traits::input_parameter< bool >::type non_negative(non_negativeSEXP); - initialize_biases_double(m_csc_r, m_csr_r, user_bias, item_bias, lambda, non_negative); - return R_NilValue; + rcpp_result_gen = Rcpp::wrap(initialize_biases_double(m_csc_r, m_csr_r, user_bias, item_bias, lambda, non_negative)); + return rcpp_result_gen; END_RCPP } // initialize_biases_float -void initialize_biases_float(const Rcpp::S4& m_csc_r, const Rcpp::S4& m_csr_r, Rcpp::S4& user_bias, Rcpp::S4& item_bias, double lambda, bool non_negative); +double initialize_biases_float(const Rcpp::S4& m_csc_r, const Rcpp::S4& m_csr_r, Rcpp::S4& user_bias, Rcpp::S4& item_bias, double lambda, bool non_negative); RcppExport SEXP _rsparse_initialize_biases_float(SEXP m_csc_rSEXP, SEXP m_csr_rSEXP, SEXP user_biasSEXP, SEXP item_biasSEXP, SEXP lambdaSEXP, SEXP non_negativeSEXP) { BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; Rcpp::RNGScope rcpp_rngScope_gen; Rcpp::traits::input_parameter< const Rcpp::S4& >::type m_csc_r(m_csc_rSEXP); Rcpp::traits::input_parameter< const Rcpp::S4& >::type m_csr_r(m_csr_rSEXP); @@ -347,8 +349,8 @@ BEGIN_RCPP Rcpp::traits::input_parameter< Rcpp::S4& >::type item_bias(item_biasSEXP); Rcpp::traits::input_parameter< double >::type lambda(lambdaSEXP); Rcpp::traits::input_parameter< bool >::type non_negative(non_negativeSEXP); - initialize_biases_float(m_csc_r, m_csr_r, user_bias, item_bias, lambda, non_negative); - return R_NilValue; + rcpp_result_gen = Rcpp::wrap(initialize_biases_float(m_csc_r, m_csr_r, user_bias, item_bias, lambda, non_negative)); + return rcpp_result_gen; END_RCPP } // deep_copy diff --git a/src/SolverAlsWrmf.cpp b/src/SolverAlsWrmf.cpp index 82e43fa..48d157f 100644 --- a/src/SolverAlsWrmf.cpp +++ b/src/SolverAlsWrmf.cpp @@ -168,7 +168,7 @@ T als_explicit(const dMappedCSC& Conf, if(lambda > 0) { if (with_biases) { auto X_no_bias = X(arma::span(1, X.n_rows - 1), arma::span::all); - auto Y_no_bias = X(arma::span(1, Y.n_rows - 1), arma::span::all); + auto Y_no_bias = Y(arma::span(1, Y.n_rows - 1), arma::span::all); loss += lambda * (accu(square(X_no_bias)) + accu(square(Y_no_bias))); } else { loss += lambda * (accu(square(X)) + accu(square(Y))); @@ -304,11 +304,20 @@ double als_explicit_float(const Rcpp::S4 &m_csc_r, } template -void initialize_biases(const dMappedCSC& ConfCSC, - const dMappedCSC& ConfCSR, - arma::Col& user_bias, - arma::Col& item_bias, - T lambda, bool non_negative) { +double initialize_biases(const dMappedCSC& ConfCSC, + const dMappedCSC& ConfCSR, + arma::Col& user_bias, + arma::Col& item_bias, + T lambda, bool non_negative) { + /* Robust mean calculation */ + double glob_mean = 0; + for (size_t ix = 0; ix < ConfCSC.nnz; ix++) + glob_mean += (ConfCSC.values[ix] - glob_mean) / (double)(ix+1); + for (size_t ix = 0; ix < ConfCSC.nnz; ix++) { + ConfCSC.values[ix] -= glob_mean; + ConfCSR.values[ix] -= glob_mean; + } + for (int iter = 0; iter < 5; iter++) { item_bias.zeros(); for (int col = 0; col < ConfCSC.n_cols; col++) { @@ -330,37 +339,38 @@ void initialize_biases(const dMappedCSC& ConfCSC, user_bias[row] = std::fmax(0., user_bias[row]); } } + return glob_mean; } // [[Rcpp::export]] -void initialize_biases_double(const Rcpp::S4 &m_csc_r, - const Rcpp::S4 &m_csr_r, - arma::Col& user_bias, - arma::Col& item_bias, - double lambda, bool non_negative) { +double initialize_biases_double(const Rcpp::S4 &m_csc_r, + const Rcpp::S4 &m_csr_r, + arma::Col& user_bias, + arma::Col& item_bias, + double lambda, bool non_negative) { const dMappedCSC ConfCSC = extract_mapped_csc(m_csc_r); const dMappedCSC ConfCSR = extract_mapped_csc(m_csr_r); - initialize_biases(ConfCSC, ConfCSR, - user_bias, item_bias, - lambda, non_negative); + return initialize_biases(ConfCSC, ConfCSR, + user_bias, item_bias, + lambda, non_negative); } // [[Rcpp::export]] -void initialize_biases_float(const Rcpp::S4 &m_csc_r, - const Rcpp::S4 &m_csr_r, - Rcpp::S4& user_bias, - Rcpp::S4& item_bias, - double lambda, bool non_negative) { +double initialize_biases_float(const Rcpp::S4 &m_csc_r, + const Rcpp::S4 &m_csr_r, + Rcpp::S4& user_bias, + Rcpp::S4& item_bias, + double lambda, bool non_negative) { const dMappedCSC ConfCSC = extract_mapped_csc(m_csc_r); const dMappedCSC ConfCSR = extract_mapped_csc(m_csr_r); arma::Col user_bias_arma = exctract_float_matrix(user_bias); arma::Col item_bias_arma = exctract_float_matrix(item_bias); - initialize_biases(ConfCSC, ConfCSR, - user_bias_arma, - item_bias_arma, - lambda, non_negative); + return initialize_biases(ConfCSC, ConfCSR, + user_bias_arma, + item_bias_arma, + lambda, non_negative); } // [[Rcpp::export]] From 96c848f8203aba68bd89db043f9ba63c9a8fea2b Mon Sep 17 00:00:00 2001 From: David Cortes Date: Sun, 29 Nov 2020 23:05:37 +0200 Subject: [PATCH 2/2] add simd just in case it doesn't vectorize on its own --- src/SolverAlsWrmf.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/src/SolverAlsWrmf.cpp b/src/SolverAlsWrmf.cpp index 48d157f..a000bce 100644 --- a/src/SolverAlsWrmf.cpp +++ b/src/SolverAlsWrmf.cpp @@ -313,6 +313,7 @@ double initialize_biases(const dMappedCSC& ConfCSC, double glob_mean = 0; for (size_t ix = 0; ix < ConfCSC.nnz; ix++) glob_mean += (ConfCSC.values[ix] - glob_mean) / (double)(ix+1); +#pragma omp simd for (size_t ix = 0; ix < ConfCSC.nnz; ix++) { ConfCSC.values[ix] -= glob_mean; ConfCSR.values[ix] -= glob_mean;