如何对R中的列进行滚动求和？

Question

roll_sum 和许多其他方法（例如 https://vandomed.github.io/moving_averages.html）仅用于对行求和。我有一个很大的矩阵，我没有足够的内存来转置它。有没有办法可以直接对列进行 roll_sum ？

例如：

library(roll)

A=matrix(rnorm(10000),100)
roll_sum(A,3)

但我想跨列执行此操作。

---

接下来，到目前为止所有的方法都是在不使用多核处理的情况下实现的。任何人都可以提供具有此功能的解决方案吗？

Answer 1

这是一个rcpp方法。

Rcpp::cppFunction("
NumericMatrix rcpp_column_roll(const NumericMatrix mat, const int n) {

  const int ncol = mat.ncol();
  const int nrow = mat.nrow();
  NumericMatrix out(nrow, ncol);
  std::fill( out.begin(), out.end(), NumericVector::get_na() ) ;

  
  for (int i = 0; i < nrow; i++) {
    NumericVector window(n);
    double roll = 0;
    int oldest_ind = 0;
    
    for (int j = 0; j < n ; j++) {
      double mat_ij = mat(i, j); 
      window(j) = mat_ij;
      roll += mat_ij;
    }
    
    out(i, n - 1) = roll;

    for (int j = n; j < ncol; j ++) {
      double mat_ij = mat(i, j); 
      
      roll += mat_ij;
      roll -= window(oldest_ind);
      
      out(i, j) = roll;
      
      window(oldest_ind) = mat_ij;
      
      if (oldest_ind == n-1) oldest_ind = 0; else oldest_ind++;
    }
  }
  return(out);
}
")

这比转置apply(A, 1L, roll::roll_sum, 3L)样本数据集的结果大约高 10 倍的内存效率和大约 50 倍的速度。

bench::mark(rcpp_column_roll(A, 3),
            t(apply(A, 1, roll::roll_sum, 3)))

## # A tibble: 2 x 13
##   expression                             min   median `itr/sec` mem_alloc
##   <bch:expr>                        <bch:tm> <bch:tm>     <dbl> <bch:byt>
## 1 rcpp_column_roll(A, 3)             134.4us  139.7us     6641.    80.7KB
## 2 t(apply(A, 1, roll::roll_sum, 3))   7.62ms   8.91ms      101.     773KB

## With an 80 MB dataset (`rnorm(1E7)`):

##   expression                          min median `itr/sec` mem_alloc
##   <bch:expr>                        <bch> <bch:>     <dbl> <bch:byt>
## 1 rcpp_column_roll(A, 3)            226ms  229ms      4.17    76.3MB
## 2 t(apply(A, 1, roll::roll_sum, 3)) 740ms  740ms      1.35   498.5MB

## 800 MB dataset (`rnorm(1E8)`):

## # A tibble: 2 x 13
##   expression                          min median `itr/sec` mem_alloc
##   <bch:expr>                        <bch> <bch:>     <dbl> <bch:byt>
## 1 rcpp_column_roll(A, 3)            3.49s  3.49s     0.286  762.94MB
## 2 t(apply(A, 1, roll::roll_sum, 3)) 9.62s  9.62s     0.104    4.84GB

内存节省似乎稳定在减少 5 倍左右，并且或多或少是结果矩阵本身的分配。

或者，我们可以更接近 R 并使用 R 循环来制作apply不需要转置的手册。

out = matrix(NA_real_, nrow(A), ncol(A))
for (i in seq_len(nrow(A))) {
  out[i, ] = roll::roll_sum(A[i, ], 3L)
}

Is 比移调常规apply. @Moody_Mudskipper 拥有最快的方法，尽管rcpp是内存效率最高的。

##rnorm(1e8); ncols = 1000;
# A tibble: 6 x 13
  expression               min median `itr/sec` mem_alloc `gc/sec` n_itr
  <bch:expr>             <bch> <bch:>     <dbl> <bch:byt>    <dbl> <int>
1 rcpp_column_roll(A, 3) 3.32s  3.32s     0.301  762.94MB    0         1
2 for_loop               6.12s  6.12s     0.163    2.98GB    0.327     1
3 dww_sappy                 7s     7s     0.143    4.86GB    0.572     1
4 matStat_Moody          1.81s  1.81s     0.552    2.24GB    0.552     1
5 roll_sum_Ronak         8.34s  8.34s     0.120    4.84GB    0.360     1
6 froll_Oliver           7.75s  7.75s     0.129    4.86GB    0.516     1

请注意，如果您的 RAM 确实不足，您可以更改 Rcpp 函数以直接修改输入，这意味着您不必分配另一个矩阵。否则，您最好在 Rcpp 中实施穆迪巧妙的解决方案，因为它会更快，并且只需要分配输出矩阵。

Answer 2

也许，您可以尝试apply按行使用矩阵：

apply(A, 1, zoo::rollsumr, 3, fill = NA)
#Or
#apply(A, 1, roll::roll_sum, 3)

However, note that this will give you output in column-order format. For example,

A <- matrix(1:10, ncol = 5)
apply(A, 1, zoo::rollsumr, 3, fill = NA)

#     [,1] [,2]
#[1,]   NA   NA
#[2,]   NA   NA
#[3,]    9   12
#[4,]   15   18
#[5,]   21   24