问题描述
我正在尝试在 Rcpp 中调用 R 函数optim().我在在R ++中从C ++内调用R的优化函数中看到了一个示例，但是对于我的用例，我无法对其进行正确的修改.基本上，目标函数取决于x和y，但我想针对b对其进行优化.
I am trying to call R function optim() in Rcpp. I saw an example in Calling R's optim function from within C++ using Rcpp, but I am unable to modify it correctly for my use case. Basically, the objective function depends on the x and y but I want to optimize it with respect to b.
这是 R 代码，可以满足我的要求:
Here is the R code that does what I want:
example_r = function(b, x, y) {
phi = rnorm(length(x))
tar_val = (x ^ 2 + y ^ 2) * b * phi
objftn_r = function(beta, x, y) {
obj_val = (x ^ 2 + y ^ 2) * beta
return(obj_val)
}
b1 = optim(b, function(beta) {
sum((objftn_r(beta, x, y) - tar_val) ^ 2)
}, method = "BFGS")$par
result = (x ^ 2 + y ^ 2) * b1
return(b1)
}
这是我尝试将其翻译为_RcppArmadillo:
Here's is my attempt to translate it to _RcppArmadillo:
#include <RcppArmadillo.h>
using namespace Rcpp;
// [[Rcpp::depends(RcppArmadillo)]]
arma::vec example_rcpp(arma::vec b, arma::vec x, arma::vec y){
arma::vec tar_val = pow(x,2)%b-pow(y,2);
return tar_val;
}
// [[Rcpp::export]]
arma::vec optim_rcpp(const arma::vec& init_val, arma::vec& x, arma::vec& y){
Rcpp::Environment stats("package:stats");
Rcpp::Function optim = stats["optim"];
Rcpp::List opt_results = optim(Rcpp::_["par"] = init_val,
Rcpp::_["fn"] = Rcpp::InternalFunction(&example_rcpp),
Rcpp::_["method"] = "BFGS");
arma::vec out = Rcpp::as<arma::vec>(opt_results[0]);
return out;
}
但是，此代码返回:
> optim_rcpp(1:3,2:4,3:5) Error in optim_rcpp(1:3, 2:4, 3:5) : not compatible with requested type
我不确定这是什么错误.
I'm not sure what the error is here.
推荐答案
在开始之前，我有几点评论:
Before we begin, I have a few remarks:
- 请显示您的所有尝试.
- 尤其要确保您的示例是可重现的最小示例
- Please show all of your attempt.
- In particular, make sure your example is a minimal reproducible example
- 在 C ++ 中使用 R 中的optim与在 C ++ 中使用基础 C ++ 代码非常不同用于nlopt中的opt().
- Using optim from R in C++ is very different than using in C++ the underlying C++ code for opt() from nlopt.
- 如果您连续快速地问了三个以上的问题，请阅读文档或与熟悉此内容的人面谈.
结果，我已经清理了您的问题...但是，将来可能不会发生这种情况.
I've cleaned up your question as a result... But, in the future, this likely will not happen.
数据生成过程似乎分两个步骤完成:首先在example_r函数外部，然后在函数内部.
The data generation process seems to be done in 2 steps: First, outside of the example_r function, and, then inside the function.
应对此进行简化，以使其在 optimization 函数之外进行.例如:
This should be simplified so that it is done outside of the optimization function. For example:
generate_data = function(n, x_mu = 0, y_mu = 1, beta = 1.5) {
x = rnorm(n, x_mu)
y = rnorm(n, y_mu)
phi = rnorm(length(x))
tar_val = (x ^ 2 + y ^ 2) * beta * phi
simulated_data = list(x = x, y = y, beta = beta, tar_val = tar_val)
return(simulated_data)
}
目标函数和 R 的optim
目标函数必须返回单个值，例如 R 中的标量.在发布的 R 代码下，实际上有两个功能被设计为按顺序充当目标功能，例如
Objective Functions and R's optim
Objective functions must return a single value, e.g. a scalar, in R. Under the posted R code, there was effectively two functions designed to act as an objective function in sequence, e.g.
objftn_r = function(beta, x, y) {
obj_val = (x ^ 2 + y ^ 2) * beta
return(obj_val)
}
b1 = optim(b, function(beta) {
sum((objftn_r(beta, x, y) - tar_val) ^ 2)
}, method = "BFGS")$par
因此，该目标函数应重写为:
This objective function should therefore be re-written as:
objftn_r = function(beta_hat, x, y, tar_val) {
# The predictions generate will be a vector
est_val = (x ^ 2 + y ^ 2) * beta_hat
# Here we apply sum of squares which changes it
# from a vector into a single "objective" value
# that optim can work with.
obj_val = sum( ( est_val - tar_val) ^ 2)
return(obj_val)
}
从那里，呼叫应对齐为:
From there, the calls should align as:
sim_data = generate_data(10, 1, 2, .3)
b1 = optim(sim_data$beta, fn = objftn_r, method = "BFGS",
x = sim_data$x, y = sim_data$y, tar_val = sim_data$tar_val)$par
RcppArmadillo目标函数
固定了 R 代码的范围和行为后，我们集中精力将其翻译为 RcppArmadillo .
RcppArmadillo Objective Functions
Having fixed the scope and behavior of the R code, let's focus on translating it into RcppArmadillo.
尤其要注意，转换后定义的异议函数将 vector 而不是 scalar 返回到optim，而 not 单个值.另外值得关注的是目标函数中缺少tar_val参数.考虑到这一点，目标函数将转换为:
In particular, notice that the objection function defined after the translation returns a vector and not a scalar into optim, which is not a single value. Also of concern is the lack of a tar_val parameter in the objective function. With this in mind, the objective function will translate to:
// changed function return type and
// the return type of first parameter
double obj_fun_rcpp(double& beta_hat,
arma::vec& x, arma::vec& y, arma::vec& tar_val){
// Changed from % to * as it is only appropriate if
// `beta_hat` is the same length as x and y.
// This is because it performs element-wise multiplication
// instead of a scalar multiplication on a vector
arma::vec est_val = (pow(x, 2) - pow(y, 2)) * beta_hat;
// Compute objective value
double obj_val = sum( pow( est_val - tar_val, 2) );
// Return a single value
return obj_val;
}
现在，在设置了目标函数的情况下，让我们针对 C ++ 中optim()的 R 调用 Rcpp .在此功能中，必须明确提供功能 .因此，x，y和tar_val必须出现在optim调用中.因此，我们将得出以下结论:
Now, with the objective function set, let's address the Rcpp call into R for optim() from C++. In this function, the parameters of thefunction must be explicitly supplied. So, x, y, and tar_val must be present in the optim call. Thus, we will end up with:
// [[Rcpp::export]]
arma::vec optim_rcpp(double& init_val,
arma::vec& x, arma::vec& y, arma::vec& tar_val){
// Extract R's optim function
Rcpp::Environment stats("package:stats");
Rcpp::Function optim = stats["optim"];
// Call the optim function from R in C++
Rcpp::List opt_results = optim(Rcpp::_["par"] = init_val,
// Make sure this function is not exported!
Rcpp::_["fn"] = Rcpp::InternalFunction(&obj_fun_rcpp),
Rcpp::_["method"] = "BFGS",
// Pass in the other parameters as everything
// is scoped environmentally
Rcpp::_["x"] = x,
Rcpp::_["y"] = y,
Rcpp::_["tar_val"] = tar_val);
// Extract out the estimated parameter values
arma::vec out = Rcpp::as<arma::vec>(opt_results[0]);
// Return estimated values
return out;
}
在一起
完整功能代码可以用test_optim.cpp编写，并通过sourceCpp()编译为:
All together
The full functioning code can be written in test_optim.cpp and compiled via sourceCpp() as:
#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
// changed function return type and
// the return type of first parameter
// DO NOT EXPORT THIS FUNCTION VIA RCPP ATTRIBUTES
double obj_fun_rcpp(double& beta_hat,
arma::vec& x, arma::vec& y, arma::vec& tar_val){
// Changed from % to * as it is only appropriate if
// `beta_hat` is the same length as x and y.
// This is because it performs element-wise multiplication
// instead of a scalar multiplication on a vector
arma::vec est_val = (pow(x, 2) - pow(y, 2)) * beta_hat;
// Compute objective value
double obj_val = sum( pow( est_val - tar_val, 2) );
// Return a single value
return obj_val;
}
// [[Rcpp::export]]
arma::vec optim_rcpp(double& init_val,
arma::vec& x, arma::vec& y, arma::vec& tar_val){
// Extract R's optim function
Rcpp::Environment stats("package:stats");
Rcpp::Function optim = stats["optim"];
// Call the optim function from R in C++
Rcpp::List opt_results = optim(Rcpp::_["par"] = init_val,
// Make sure this function is not exported!
Rcpp::_["fn"] = Rcpp::InternalFunction(&obj_fun_rcpp),
Rcpp::_["method"] = "BFGS",
// Pass in the other parameters as everything
// is scoped environmentally
Rcpp::_["x"] = x,
Rcpp::_["y"] = y,
Rcpp::_["tar_val"] = tar_val);
// Extract out the estimated parameter values
arma::vec out = Rcpp::as<arma::vec>(opt_results[0]);
// Return estimated values
return out;
}
测试用例
# Setup some values beta = 2 x = 2:4 y = 3:5 # Set a seed for reproducibility set.seed(111) phi = rnorm(length(x)) tar_val = (x ^ 2 + y ^ 2) * beta * phi optim_rcpp(beta, x, y, tar_val) # [,1] # [1,] 2.033273
注意:如果要避免返回大小为1 x1的矩阵，请使用double作为optim_rcpp的返回参数，并将Rcpp::as<arma::vec>切换为Rcpp::as<double>
Note: If you would like to avoid a matrix of size 1 x1 from being returned please use double as the return parameter of optim_rcpp and switch Rcpp::as<arma::vec> to Rcpp::as<double>
这篇关于使用Rcpp在C ++中的R中应用优化函数的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持！