Merge pull request #296 from xoopR/fix_weightvec_bottleneck

fix massive bottleneck
xoopR · Nov 11, 2023 · a656b7b · a656b7b
2 parents 255a666 + ddef09a
commit a656b7b
Show file tree

Hide file tree

Showing 7 changed files with 45 additions and 51 deletions.
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: distr6
 Title: The Complete R6 Probability Distributions Interface
-Version: 1.8.3
+Version: 1.8.4
 Authors@R:
     c(person(given = "Raphael",
              family = "Sonabend",

diff --git a/NEWS.md b/NEWS.md
@@ -1,3 +1,7 @@
+# distr6 1.8.4
+
+* Fix massive bottleneck in Matdist and Arrdist
+
 # distr6 1.8.3
 
 * Add `decorators` argument to `c.Matdist` and `c.Arrdist`

diff --git a/R/SDistribution_Arrdist.R b/R/SDistribution_Arrdist.R
@@ -310,30 +310,28 @@ Arrdist <- R6Class("Arrdist",
     .pdf = function(x, log = FALSE) {
       "pdf, data, wc" %=% gprm(self, c("pdf", "x", "which.curve"))
       mat <- .extCurve(pdf, wc)
-      out <- t(C_Vec_WeightedDiscretePdf(
-        x, matrix(data, ncol(mat), private$.ndists), t(mat)))
+      out <- t(C_Vec_WeightedDiscretePdf(x, data, t(mat)))
       if (log) {
         out <- log(out)
       }
       colnames(out) <- x
       t(out)
     },
 
-    .cdf = function(x, lower.tail = TRUE, log.p = FALSE) { # FIXME
+    .cdf = function(x, lower.tail = TRUE, log.p = FALSE) {
       "cdf, data, wc" %=% gprm(self, c("cdf", "x", "which.curve"))
       mat <- .extCurve(cdf, wc)
-      out <- t(C_Vec_WeightedDiscreteCdf(
-        x, matrix(data, ncol(mat), nrow(mat)), t(mat), lower.tail, log.p
-      ))
+      out <- t(C_Vec_WeightedDiscreteCdf(x, data, t(mat), lower.tail,
+        log.p))
       colnames(out) <- x
       t(out)
     },
 
     .quantile = function(p, lower.tail = TRUE, log.p = FALSE) {
       "*" %=% gprm(self, c("cdf", "x", "which.curve"))
       mat <- .extCurve(cdf, which.curve)
-      out <- t(C_Vec_WeightedDiscreteQuantile(p,
-        matrix(x, ncol(mat), nrow(mat)), t(mat), lower.tail, log.p))
+      out <- t(C_Vec_WeightedDiscreteQuantile(p, x, t(mat), lower.tail,
+        log.p))
       colnames(out) <- NULL
       t(out)
     },

diff --git a/R/SDistribution_Matdist.R b/R/SDistribution_Matdist.R
@@ -274,31 +274,26 @@ Matdist <- R6Class("Matdist",
     # dpqr
     .pdf = function(x, log = FALSE) {
       "pdf, data" %=% gprm(self, c("pdf", "x"))
-      out <- t(C_Vec_WeightedDiscretePdf(
-        x, matrix(data, ncol(pdf), private$.ndists),
-        t(pdf)
-      ))
+      out <- t(C_Vec_WeightedDiscretePdf(x, data, t(pdf)))
       if (log) {
         out <- log(out)
       }
       colnames(out) <- x
       t(out)
     },
 
-    .cdf = function(x, lower.tail = TRUE, log.p = FALSE) { # FIXME
+    .cdf = function(x, lower.tail = TRUE, log.p = FALSE) {
       "cdf, data" %=% gprm(self, c("cdf", "x"))
-      out <- t(C_Vec_WeightedDiscreteCdf(
-        x, matrix(data, ncol(cdf), nrow(cdf)),
-        t(cdf), lower.tail, log.p
-      ))
+      out <- t(C_Vec_WeightedDiscreteCdf(x, data, t(cdf), lower.tail,
+        log.p))
       colnames(out) <- x
       t(out)
     },
 
     .quantile = function(p, lower.tail = TRUE, log.p = FALSE) {
       "*" %=% gprm(self, c("cdf", "x"))
-      out <- t(C_Vec_WeightedDiscreteQuantile(p, matrix(x, ncol(cdf), nrow(cdf)),
-          t(cdf), lower.tail, log.p))
+      out <- t(C_Vec_WeightedDiscreteQuantile(p, x, t(cdf), lower.tail,
+        log.p))
       colnames(out) <- NULL
       t(out)
     },

diff --git a/R/SDistribution_WeightedDiscrete.R b/R/SDistribution_WeightedDiscrete.R
@@ -375,8 +375,7 @@ WeightedDiscrete <- R6Class("WeightedDiscrete",
           data[[i]] <- data[[i]][seq.int(lng)]
         }
         pdf <- matrix(unlist(pdf), nrow = length(data[[1]]), ncol = length(data))
-        data <- matrix(unlist(data), ncol = ncol(pdf))
-        out <- C_Vec_WeightedDiscretePdf(x, data, pdf)
+        out <- C_Vec_WeightedDiscretePdf(x, unlist(data), pdf)
         if (log) {
           out <- log(out)
         }
@@ -396,8 +395,7 @@ WeightedDiscrete <- R6Class("WeightedDiscrete",
           data[[i]] <- data[[i]][seq.int(lng)]
         }
         cdf <- matrix(unlist(cdf), nrow = length(data[[1]]), ncol = length(data))
-        data <- matrix(unlist(data), ncol = ncol(cdf))
-        C_Vec_WeightedDiscreteCdf(x, data, cdf, lower.tail, log.p)
+        C_Vec_WeightedDiscreteCdf(x, unlist(data), cdf, lower.tail, log.p)
       } else {
         .wd_cdf(x, data, cdf, lower.tail, log.p)
       }
@@ -414,8 +412,7 @@ WeightedDiscrete <- R6Class("WeightedDiscrete",
           data[[i]] <- data[[i]][seq.int(lng)]
         }
         cdf <- matrix(unlist(cdf), nrow = length(data[[1]]), ncol = length(data))
-        data <- matrix(unlist(data), ncol = ncol(cdf))
-        C_Vec_WeightedDiscreteQuantile(p, data, cdf, lower.tail, log.p)
+        C_Vec_WeightedDiscreteQuantile(p, unlist(data), cdf, lower.tail, log.p)
       } else {
         C_WeightedDiscreteQuantile(p, data, cdf, lower.tail, log.p)
       }

diff --git a/src/Distributions.cpp b/src/Distributions.cpp
@@ -399,23 +399,23 @@ NumericMatrix C_ShiftedLoglogisticQuantile(NumericVector x, NumericVector locati
 
 
 // [[Rcpp::export]]
-NumericMatrix C_Vec_WeightedDiscretePdf(NumericVector x, NumericMatrix data,
+NumericMatrix C_Vec_WeightedDiscretePdf(NumericVector x, NumericVector data,
                                         NumericMatrix pdf) {
 
-  int nc = data.ncol();
-  int nr = data.nrow();
+  int nc = pdf.ncol();
+  int nr = pdf.nrow();
   int n = x.length();
 
   NumericMatrix mat(n, nc);
 
-  // i - distribution
-  // j - data samples
-  // k - evaluates
+  // i - data samples
+  // j - data points
+  // k - new data points
 
   for (int i = 0; i < nc; i++) {
     for (int k = 0; k < n; k++) {
       for (int j = 0; j < nr; j++) {
-        if (data(j, i) == x[k]) {
+        if (data[j] == x[k]) {
           mat(k, i) = pdf(j, i);
           break;
         }
@@ -467,29 +467,29 @@ NumericVector C_WeightedDiscreteCdf(NumericVector x, NumericVector data, Numeric
 }
 
 // [[Rcpp::export]]
-NumericMatrix C_Vec_WeightedDiscreteCdf(NumericVector x, NumericMatrix data, NumericMatrix cdf,
+NumericMatrix C_Vec_WeightedDiscreteCdf(NumericVector x, NumericVector data, NumericMatrix cdf,
                                         bool lower, bool logp) {
 
-  int nc = data.ncol();
-  int nr = data.nrow();
+  int nc = cdf.ncol();
+  int nr = cdf.nrow();
   int n = x.length();
 
   NumericMatrix mat(n, nc);
 
-  // i - distribution
-  // j - data samples
-  // k - evaluates
+  // i - data samples
+  // j - data points
+  // k - new data points
 
   for (int i = 0; i < nc; i++) {
     for (int k = 0; k < n; k++) {
       for (int j = 0; j < nr; j++) {
-        if (j == 0 && x[k] < data(0, i)) {
+        if (j == 0 && x[k] < data[0]) {
           mat(k, i) = 0;
           break;
         } else if (j == nr - 1) {
           mat(k, i) = cdf(j, i);
           break;
-        } else if (x[k] >= data(j, i) && x[k] < data(j + 1, i)) {
+        } else if (x[k] >= data[j] && x[k] < data[j + 1]) {
           mat(k, i) = cdf(j, i);
           break;
         }
@@ -540,11 +540,11 @@ NumericVector C_WeightedDiscreteQuantile(NumericVector x, NumericVector data, Nu
 }
 
 // [[Rcpp::export]]
-NumericMatrix C_Vec_WeightedDiscreteQuantile(NumericVector x, NumericMatrix data, NumericMatrix cdf,
+NumericMatrix C_Vec_WeightedDiscreteQuantile(NumericVector x, NumericVector data, NumericMatrix cdf,
                                              bool lower, bool logp) {
 
-  int nc = data.ncol();
-  int nr = data.nrow();
+  int nc = cdf.ncol();
+  int nr = cdf.nrow();
   int n = x.length();
 
   NumericMatrix mat(n, nc);
@@ -567,7 +567,7 @@ NumericMatrix C_Vec_WeightedDiscreteQuantile(NumericVector x, NumericMatrix data
         }
 
         if (y <= cdf(j, i)) {
-          mat(k, i) = data(j, i);
+          mat(k, i) = data[j];
           break;
         }
       }

diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp
@@ -179,13 +179,13 @@ BEGIN_RCPP
 END_RCPP
 }
 // C_Vec_WeightedDiscretePdf
-NumericMatrix C_Vec_WeightedDiscretePdf(NumericVector x, NumericMatrix data, NumericMatrix pdf);
+NumericMatrix C_Vec_WeightedDiscretePdf(NumericVector x, NumericVector data, NumericMatrix pdf);
 RcppExport SEXP _distr6_C_Vec_WeightedDiscretePdf(SEXP xSEXP, SEXP dataSEXP, SEXP pdfSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
     Rcpp::traits::input_parameter< NumericVector >::type x(xSEXP);
-    Rcpp::traits::input_parameter< NumericMatrix >::type data(dataSEXP);
+    Rcpp::traits::input_parameter< NumericVector >::type data(dataSEXP);
     Rcpp::traits::input_parameter< NumericMatrix >::type pdf(pdfSEXP);
     rcpp_result_gen = Rcpp::wrap(C_Vec_WeightedDiscretePdf(x, data, pdf));
     return rcpp_result_gen;
@@ -207,13 +207,13 @@ BEGIN_RCPP
 END_RCPP
 }
 // C_Vec_WeightedDiscreteCdf
-NumericMatrix C_Vec_WeightedDiscreteCdf(NumericVector x, NumericMatrix data, NumericMatrix cdf, bool lower, bool logp);
+NumericMatrix C_Vec_WeightedDiscreteCdf(NumericVector x, NumericVector data, NumericMatrix cdf, bool lower, bool logp);
 RcppExport SEXP _distr6_C_Vec_WeightedDiscreteCdf(SEXP xSEXP, SEXP dataSEXP, SEXP cdfSEXP, SEXP lowerSEXP, SEXP logpSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
     Rcpp::traits::input_parameter< NumericVector >::type x(xSEXP);
-    Rcpp::traits::input_parameter< NumericMatrix >::type data(dataSEXP);
+    Rcpp::traits::input_parameter< NumericVector >::type data(dataSEXP);
     Rcpp::traits::input_parameter< NumericMatrix >::type cdf(cdfSEXP);
     Rcpp::traits::input_parameter< bool >::type lower(lowerSEXP);
     Rcpp::traits::input_parameter< bool >::type logp(logpSEXP);
@@ -237,13 +237,13 @@ BEGIN_RCPP
 END_RCPP
 }
 // C_Vec_WeightedDiscreteQuantile
-NumericMatrix C_Vec_WeightedDiscreteQuantile(NumericVector x, NumericMatrix data, NumericMatrix cdf, bool lower, bool logp);
+NumericMatrix C_Vec_WeightedDiscreteQuantile(NumericVector x, NumericVector data, NumericMatrix cdf, bool lower, bool logp);
 RcppExport SEXP _distr6_C_Vec_WeightedDiscreteQuantile(SEXP xSEXP, SEXP dataSEXP, SEXP cdfSEXP, SEXP lowerSEXP, SEXP logpSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
     Rcpp::traits::input_parameter< NumericVector >::type x(xSEXP);
-    Rcpp::traits::input_parameter< NumericMatrix >::type data(dataSEXP);
+    Rcpp::traits::input_parameter< NumericVector >::type data(dataSEXP);
     Rcpp::traits::input_parameter< NumericMatrix >::type cdf(cdfSEXP);
     Rcpp::traits::input_parameter< bool >::type lower(lowerSEXP);
     Rcpp::traits::input_parameter< bool >::type logp(logpSEXP);