从 rcpp 返回 R 函数

Question

Rcpp 中有没有办法返回一个 R 函数，其中包含一些仅在第一次函数调用时计算的预先计算的值？请考虑以下 R 代码：

1: func_generator<-function(X) {
2:  X_tot<-sum(X)
3:  function(b_vec) { (X_tot*b_vec) }
4: }
5: myfunc<-func_generator(c(3,4,5))
6: myfunc(1:2)
7: myfunc(5:6)
8: myfunc2<-func_generator(c(10,11,12,13))
...

这可以在 Rcpp 中编程吗？在实践中，假设用更密集的计算来代替第 2 行。

为了添加上下文，给定向量 X 和标量 b，存在一些似然函数 f(b|X(，对于一些足够的统计量s(X(，它可以重新表示为f(b，s(X((，该统计量仅是X的函数，并且涉及一些计算。这是在计算密集型计算机实验中，有许多向量X(许多可能性(，以及许多对f(bvec|X(对于每个可能性，所以我宁愿计算一次s(X((对于每个可能性(并以某种方式保存它，而不是多次重新计算它。我从简单地编程 f(bvec，X( 开始计算 f(b|X( 在点 bvec=(b_1,...,b_n(，但这有额外的开销，因为我多次调用此函数并在每次运行时计算 s(X(。我想只计算一次 s(X(。

任何在 Rcpp 中有效完成此任务的建议都将不胜感激(无论是通过返回函数;还是通过以其他方式存储中间计算(。

Answer 1

存储中间结果的一种简单方法是在函数级别使用静态变量：

// [[Rcpp::plugins(cpp11)]]
#include <thread>
#include <chrono>
#include <Rcpp.h>
// [[Rcpp::export]]
Rcpp::NumericVector foo(Rcpp::NumericVector X, Rcpp::NumericVector b, bool useCache = true) {
static double cache;
static bool initialized{false};
if (!(useCache && initialized)) {
// sleep to simulate actual work
std::this_thread::sleep_for (std::chrono::seconds(1));
cache = Rcpp::sum(X);
initialized = true;
}
return cache * b;
}
/*** R
X <- 1:10
b <- 10:20
system.time(r1 <- foo(X, b))
system.time(r2 <- foo(X, b))
all.equal(r1, r2)
system.time(r3 <- foo(X, b, FALSE))
all.equal(r1, r3)
*/

输出：

> system.time(r1 <- foo(X, b))
user  system elapsed 
0       0       1 
> system.time(r2 <- foo(X, b))
user  system elapsed 
0.002   0.000   0.002 
> all.equal(r1, r2)
[1] TRUE
> system.time(r3 <- foo(X, b, FALSE))
user  system elapsed 
0       0       1 
> all.equal(r1, r3)
[1] TRUE

在第二个函数调用中使用缓存时，几乎是即时计算结果的。

如果可以在不同X的循环中循环不同的b，则此方法非常有效。如果此限制对您不起作用，则可以在 R 级别使用memoise包来有效地存储任意输入的昂贵函数的输出：

// [[Rcpp::plugins(cpp11)]]
#include <thread>
#include <chrono>
#include <Rcpp.h>
// [[Rcpp::export]]
Rcpp::NumericVector foo(double total, Rcpp::NumericVector b) {
return total * b;
}
// [[Rcpp::export]]
double bar(Rcpp::NumericVector X) {
// sleep to simulate actual work
std::this_thread::sleep_for (std::chrono::seconds(1));
return Rcpp::sum(X);
}

/*** R
X1 <- 1:10
b1 <- 10:20
X2 <- 10:1
b2 <- 20:10
library(memoise)
bar2 <- memoise(bar)
system.time(r11 <- foo(bar2(X1), b1))
system.time(r21 <- foo(bar2(X2), b2))
system.time(r12 <- foo(bar2(X1), b1))
system.time(r22 <- foo(bar2(X2), b2))
all.equal(r11, r12)
all.equal(r21, r22)
*/

输出：

> system.time(r11 <- foo(bar2(X1), b1))
user  system elapsed 
0.001   0.000   1.001 
> system.time(r21 <- foo(bar2(X2), b2))
user  system elapsed 
0.033   0.000   1.033 
> system.time(r12 <- foo(bar2(X1), b1))
user  system elapsed 
0       0       0 
> system.time(r22 <- foo(bar2(X2), b2))
user  system elapsed 
0       0       0 
> all.equal(r11, r12)
[1] TRUE
> all.equal(r21, r22)
[1] TRUE

作为替代方法，您也可以将这两个函数用作函数生成器的构建块：

func_generator <- function(X) {
X_tot <- bar(X)
function(b_vec) { foo(X_tot, b_vec) }
}
myfunc <- func_generator(c(3,4,5))
myfunc2 <- func_generator(c(10,11,12,13))
myfunc(1:2)
myfunc(5:6)
myfunc2(1:2)
myfunc2(5:6)

因此，将数字昂贵的工作保持在C++，但要保持简单。然后可以使用 R 添加功能方面。