替换R中的for循环以加速代码

Question

我使用R来调查返回如何影响个人的养老金账户。为了做到这一点，我计算了退休金账户从25岁到退休到70岁退休的1000种不同的回报情况。我使用变量费用（e），月度存款（m），百分比（r）回报，账户余额（y）和欧元回报（x）。它们都存储在尺寸为46x1000的数据框中。

我成功设法使用for循环来计算它。然而，这是非常缓慢的，因为我正在做很多这些我想知道如果有人有一个想法加快代码。我已经尝试过应用函数和矢量化，但无法使其工作。我的问题是，在计算第i + 1年的数字之前，我必须计算第i年的数字。我已经在网上搜索了一个解决方案，但是很难找到适用于我的特定问题的答案。我要指出，我还是很新的R.

一个已经使用编写的代码IM的简化版本：

for (i in 3:46) { 
x[i-1,]<-(y[i-1,]+m[i-1,]*6-0.5*e[i-1,])*r[i-1,] 
y[i,]<-y[i-1,]+x[i-1,]-e[i-1,]+m[i-1,]*12 
} 

我希望有人能够帮助，并在此先感谢。

问候 拉斯穆斯

Answer 1

你的过程在我看来就像它需要循环，因为每次迭代依赖于面前的一个。正如@Gregor de Cillia在评论中提到的那样，你可以用C++来提高速度。

首先，设置一些数据。

set.seed(1) 
e <- matrix(data = rnorm(n = 46000, mean = 1000, sd = 200), 
         nrow = 46, 
         ncol = 1000) 
m <- matrix(data = rnorm(n = 46000, mean = 2000, sd = 200), 
         nrow = 46, 
         ncol = 1000) 
r <- matrix(data = rnorm(n = 46000, mean = 4, sd = 0.5), 
         nrow = 46, 
         ncol = 1000) 
x <- matrix(data = NA_real_, nrow = 45, ncol = 1000) 
y <- matrix(data = NA_real_, nrow = 46, ncol = 1000) 
y[1,] <- rnorm(n = 1000, 10000, 1000) 

然后在Rcpp文件中定义一个C++函数。此方法返回的两个矩阵x和y列表项的列表：

List pension(NumericMatrix e, 
       NumericMatrix m, 
       NumericMatrix r, 
       NumericVector yfirstrow) { 

    int ncols = e.cols(); 
    int nrows = e.rows(); 

    NumericMatrix x(nrows - 1, ncols); 
    NumericMatrix y(nrows, ncols); 

    y(0, _) = yfirstrow; 

    for(int i = 1; i < nrows; i++) { 
     x(i-1, _) = (y(i-1, _) + m(i-1, _) * 6 - 0.5 * e(i-1, _)) * r(i-1, _); 
     y(i, _) = y(i-1, _) + x(i-1, _) - e(i-1, _) + m(i-1, _)* 12; 
    }; 

    List ret; 
    ret["x"] = x; 
    ret["y"] = y; 

    return ret; 

} 

比较对速度的两种方法。

microbenchmark::microbenchmark(
    R = { 
     for (i in 2:46) { 
      x[i-1,] <- unlist((y[i-1,] + m[i-1,]*6 - 0.5*e[i-1,]) * r[i-1,]) 
      y[i,]<- unlist(y[i-1,]+x[i-1,]-e[i-1,]+m[i-1,]*12) 
     } 
    }, 
    cpp = { 
     cppList <- pension(e, m, r, y[1,]) 
    }, 
    times = 100 
) 

确保输出匹配：

> identical(x, cppList$x) 
[1] TRUE 
> identical(y, cppList$y) 
[1] TRUE 

速度测试结果：

Unit: microseconds 
expr  min  lq  mean median  uq  max neval 
    R 3309.962 3986.569 6961.838 5244.479 6219.215 96576.592 100 
    cpp 879.713 992.229 1266.014 1124.345 1273.691 3041.966 100 

所以Rcpp解决方案是围绕更快这里5倍，但说实话，在R循环你所做的对于你正在使用的数据集来说并不是太简单（只有45次迭代，R循环的开销并不是太大的障碍）。如果你真的需要这个速度，C++可以提供帮助。