R 中的递归均值与 Rcpp

我正在尝试使用简单的递归实现来计算变量的平均值：

m <- 0 # initialize mean

for(irep in 0:999){

# new data point
new_data <- rnorm(1,2,1)

# recursive formula for sample mean
m = (irep/(irep+1)) * m + (1/(irep+1)) * new_data

}

Here, m将很快收敛到 2，它对应于我们生成新数据点的正态分布的平均值。在 Rcpp 中实现类似的东西：

#include <RcppArmadillo.h>

// [[Rcpp::depends(RcppArmadillo)]]

using namespace Rcpp;

// [[Rcpp::export]]
double  rec_mean(int sample_size){

  double m = 0; //initialize mean

  for(int irep = 0; irep < sample_size; irep++){

      // new data
      double new_data = R::rnorm(2,1);

      // mean recursive update
      m = ((irep)/(irep+1)) * m + (1/(irep+1)) * new_data;

  }

  return m;

}

此代码没有显示预期的行为。相反，它返回初始值。有人能告诉我从 R 到 Rcpp 的翻译错误在哪里吗？

在这一行中：

m = ((irep)/(irep+1)) * m + (1/(irep+1)) * new_data;

你正在划分一个int由另一个int两次。在 C++ 中，整数除法返回另一个整数，并丢弃余数。

为了得到你想要的，强制以浮点数进行除法：

m = ((irep)/(irep+1.0)) * m + (1.0/(irep+1.0)) * new_data;