将R转换为Rcpp:迭代器

马西奥(Marcio Augusto Diniz)

我正在将一个函数从R转换为Rcpp,并且已经苦苦挣扎了一段时间。我看过这样的一些材料一个,但我在RCPP一个初学者,我一直没能找出我目前做错了。

我在Rcpp中的最佳成绩如下,

#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
Rcpp::NumericVector aux_pred_C(IntegerVector m, 
                               NumericVector a, 
                               NumericVector b, 
                               int n_group, 
                               IntegerVector group_index,
                               NumericMatrix theta_sample, 
                               IntegerMatrix prior_parm, 
                               int n_mcmc){
  
  NumericVector p(n_mcmc);
  group_index = group_index - 1;

  for (int j = 0; j < n_mcmc; ++j) {
    NumericMatrix parm(2, n_group);
    NumericMatrix sample(n_mcmc, n_group);
    
    for(IntegerVector::iterator i = group_index.begin(); i != group_index.end(); ++i) {
      NumericVector temp = Rcpp::rbinom(m(i) , 1, theta_sample(j, i));
      
      parm(0,i) = prior_parm(0,i) + a(i) + sum(temp);
      parm(1,i) = prior_parm(1,i) + b(i) + m(i) - sum(temp);
      sample(_,i) = Rcpp::rbeta(n_mcmc, parm(0,i), parm(1,i));
    }

    int i1 = group_index[0];
    int i2 = group_index[1];
    LogicalVector V = sample(_,i2) > sample(_,i1);
    IntegerVector res = ifelse(V, 1, 0);
    p(j) = mean(res);
  }
  
  return(p);
}

在第22、24、25和26行重复我的错误消息:

sampler_predictive_distribution.cpp:22:44: error: no match for call to '(Rcpp::IntegerVector {aka Rcpp::Vector<13>}) (Rcpp::traits::storage_type<13>::type*&)'

invalid conversion from 'Rcpp::Vector<13>::iterator' {aka 'int*'} to 'size_t' {aka 'long long unsigned int'} [-fpermissive] NumericVector temp = Rcpp::rbinom(m(i) , 1, theta_sample(j, i));

我的猜测是我没有正确使用迭代器。我遵循了这个例子

最后,我在R中的可复制示例,

n.group <- 4
group.index <- c(1, 3)
m <- rep(NA, n.group)
m[group.index[1]] <- 50
m[group.index[2]] <- 50
a <- rep(NA, n.group)
a[group.index[1]] <- 20
a[group.index[2]] <- 20
b <- rep(NA, n.group)
b[group.index[1]] <- 30
b[group.index[2]] <- 30

n.mcmc <- 100
theta.sample.e1 <- matrix(NA, nrow = n.mcmc, ncol = n.group)
theta.sample.e1[, group.index[1]] <- rbeta(n.mcmc, a[group.index[2]], b[group.index[1]])
theta.sample.e1[, group.index[2]] <- rbeta(n.mcmc, a[group.index[2]], b[group.index[2]])
prior.parm.e1 <- matrix(1, ncol = 4, nrow = 2)

aux_pred <- function(m, a, b, n.group, group.index, theta.sample, prior.parm, n.mcmc){
  
  p <- rep(NA, n.mcmc)
  
  for (j in 1:n.mcmc){
    parm <- matrix(NA, ncol = n.group, nrow = 2)
    sample <- matrix(NA, ncol = n.group, nrow = n.mcmc)
    
    for (i in group.index){
      temp <- rbinom(m[i], size = 1, prob =  theta.sample[j, i])
      
      parm[1, i] <- prior.parm[1, i] + a[i] + sum(temp)
      parm[2, i] <- prior.parm[2, i] + b[i] + length(temp) - sum(temp)
      sample[, i] <-  rbeta(n.mcmc, parm[1, i], parm[2, i])
    }
    p[j] <- mean(sample[, group.index[2]] - sample[, group.index[1]] > 0)
  }
  
  return(p)
}

aux_pred(m, a, b, n.group, group.index, theta.sample.e1, prior.parm.e1, n.mcmc)

我想念什么?

约瑟夫·伍德

您的问题源于不使用星号取消对迭代器的引用。您将在此处的示例中注意到第28章迭代器|。每个人的Rcpp,当使用迭代器在向量的元素上循环时,为了获取特定迭代器指向的值,我们必须用星号取消引用迭代器。因此,您的代码如下所示:

for(IntegerVector::iterator i = group_index.begin(); i != group_index.end(); ++i) {
    NumericVector temp = Rcpp::rbinom(m(*i) , 1, theta_sample(j, *i));
    
    parm(0,*i) = prior_parm(0,*i) + a(*i) + sum(temp);
    parm(1,*i) = prior_parm(1,*i) + b(*i) + m(*i) - sum(temp);
    sample(_,*i) = Rcpp::rbeta(n_mcmc, parm(0,*i), parm(1,*i));
}

这会很快变得混乱。对我而言,最好使用范围基本循环

for(auto i: group_index) {
    NumericVector temp = Rcpp::rbinom(m(i) , 1, theta_sample(j, i));
    
    parm(0,i) = prior_parm(0,i) + a(i) + sum(temp);
    parm(1,i) = prior_parm(1,i) + b(i) + m(i) - sum(temp);
    sample(_,i) = Rcpp::rbeta(n_mcmc, parm(0,i), parm(1,i));
}

有一个警告,范围基于循环。也就是说,它只是因为可用的C++11,但是如果你使用的最新版本的不应该是一个令人关切的是大R这使得C++11默认。另外,您始终可以放入// [[Rcpp::plugins(cpp11)]]源文件(请参阅将C ++ 11与Rcpp一起使用的第一步)。

另外,您正在使用(检查向量边界的访问器,因此效率不如[后者更危险,并给开发人员带来责任,以确保您不访问自己不拥有的内存,但是在实践中,保证这一承诺并不困难。参见第8章第2节。rcpp给大家

本文收集自互联网,转载请注明来源。

如有侵权,请联系 [email protected] 删除。

编辑于
0

我来说两句

0 条评论
登录 后参与评论

相关文章