如何动态计算 data.table 中的行值，避免 FOR 循环答案

【问题标题】：how to dynamically calculate row values in data.table, avoiding FOR loop如何动态计算 data.table 中的行值，避免 FOR 循环
【发布时间】：2019-10-11 11:26:38
【问题描述】：

考虑以下情况：

当期 t 的广告存量取决于 t-1 的广告存量，经过某种衰减调整，以及当期的额外广告联系
当前期间 t 的采购决策取决于期间的消耗水平和 t-1 中的库存水平

将其编程为 FOR 循环非常容易。但是在运行大数据集时这会变得非常慢（例如，假设 100.000 名消费者有 2 年的每日购买记录，即有更多 1 MM 行条目）。

如何编程，例如数据表？我知道如何在 data.table 中为给定的数据集使用 shift 函数。但是这个用例中的动态计算带来了一些挑战。

我使用带有 data.frame 的 FOR 循环对此进行了编程。但是大量行的性能非常非常慢

示例背后的概念： x = 考虑为例如广告联系或按时间段的购买量等。

y = 将其视为动态变量，例如取决于 x 和 y 衰减的广告存量

decay = 来自任何类型的函数，它根据之前的 y 值和时间（通常是时间段）计算 y，根据 x - vales | 动态变化事件；这里简化为随机函数

n <- 100
DF <- data.frame(x = c(1,rep(0,n-1)), y = c(1, rep(0,n-1)), decay = c(1, rep(0,n-1)), index=rep(0,n))

set.seed(10)
for(i in 2:n){
  DF$x[i] <- sample(x=c(0:2), replace = T, size = 1, prob = c(0.8, 0.15, 0.05))
  if(DF$x[i] > 0){ DF$index[i] <- 0} else { DF$index[i] <- (DF$index[i-1] + 1) }
  DF$decay[i] <- round((DF$index[i] + 1)^-0.1, 2)
  DF$y[i] <- round((DF$x[i] + DF$y[i-1]) * DF$decay[i],2)}

plot(DF$y, type="o")

【问题讨论】：

你所做的基本上是创建一个具有一定随机影响的自回归时间序列，对吗？我认为，只需研究一个简单的自回归时间序列的创建，例如 AR(1) 过程。
我会使用/需要这种类型算法的情况并不是真正的自回归时间序列，而是用于（a）具有协变量（=多个自变量）的时间序列或时间序列建模问题的模拟ARIMA 过程不起作用（例如，危险预测）。运行模拟通常要求自变量的值取决于以前的时期。上面的程序只是超级简化，目的是给出一个代码示例。
或者在马尔可夫过程这样的情况下，过程的每个状态都依赖于前一个状态，因此需要顺序生成状态。在这里我可以想到 FOR 循环。

标签： r data.table

【解决方案1】：

矢量是你的朋友。您当前的循环每次从data.frame 中提取一个向量，这会变得很昂贵。相反，您应该：

一次性将您的 x 创建为矢量
对其他变量使用向量y, index, and decay。
将其声明为函数 - 编译器将提高性能

base_loop <- function(x) {
  y <- vector('numeric', n); y[1] <- 1
  decay <- vector('numeric', n); decay[1] <- 1
  index <- vector('integer', n)
  for(i in 2:n){
    if(x[i] > 0){index[i] <- 0} else {index[i] <- (index[i-1] + 1) }
    decay[i] <- (index[i] + 1)^-0.1
    y[i] <- (x[i] + y[i-1]) * decay[i]
  }
  data.frame(x, y, decay, index)
}

set.seed(10)
n = 1E2
x <- c(1,sample(c(0:2), replace = T, size = n-1, prob = c(0.8, 0.15, 0.05)))

DF <- base_loop(x)

这也可以很容易地转换为rcpp 循环：

// [[Rcpp::export]]
DataFrame decay_func(NumericVector x) {
  IntegerVector ind = x.size();
  NumericVector decay = x.size(); decay[1] = 1;
  NumericVector y = x.size(); y[1] = 1;



  for (int i = 0; i < x.size(); i++){
    if (x[i] > 0) {
      ind[i] = 0;
    } else {
      ind[i] = ind[i-1] + 1;
    }
    decay[i] = pow(ind[i] + 1,-0.1);
    y[i] = (x[i] + y[i-1]) * decay[i];
  }

  return DataFrame::create(Named("x") = x,
                           Named("y") = y,
                           Named("decay") = decay,
                           Named("index") = ind);
}

性能

# n = 100 
# A tibble: 4 x 13
  expression          min  median `itr/sec` mem_alloc
  <bch:expr>      <bch:t> <bch:t>     <dbl> <bch:byt>
1 OP_loop          21.5ms  21.7ms      43.9  702.62KB
2 vector_loop      11.4ms  11.8ms      82.4  100.33KB
3 compiled_vector 473.2us 483.2us    2028.     9.61KB
4 rcpp_func       412.6us 423.3us    2321.    11.27KB

# n = 10,000
# A tibble: 4 x 13
  expression          min  median `itr/sec` mem_alloc
  <bch:expr>      <bch:t> <bch:t>     <dbl> <bch:byt>
1 OP_loop           1.23s   1.23s     0.816    3.01GB
2 vector_loop     16.73ms 17.11ms    56.2    525.72KB
3 compiled_vector   5.8ms  5.88ms   167.        435KB
4 rcpp_func        1.52ms  1.55ms   606.     359.32KB

# n= 1,000,000
# A tibble: 3 x 13
  expression        min median `itr/sec` mem_alloc
  <bch:expr>      <bch> <bch:>     <dbl> <bch:byt>
1 vector_loop     563ms  563ms      1.78    42.1MB
2 compiled_vector 556ms  556ms      1.80      42MB
3 rcpp_func       115ms  120ms      7.56    34.3MB

参考代码：

library(Rcpp)

set.seed(10)
n = 1E6
x <- c(1,sample(c(0:2), replace = T, size = n-1, prob = c(0.8, 0.15, 0.05)))

bench::mark(
  # OP_loop = {
  #   DF <- data.frame(x = c(1,rep(0,n-1)), y = c(1, rep(0,n-1)), decay = c(1, rep(0,n-1)), index=rep(0,n))
  #   
  #   set.seed(10)
  #   for(i in 2:n){
  #     DF$x[i] <- sample(x=c(0:2), replace = T, size = 1, prob = c(0.8, 0.15, 0.05))
  #     if(DF$x[i] > 0){ DF$index[i] <- 0} else { DF$index[i] <- (DF$index[i-1] + 1) }
  #     DF$decay[i] <- (DF$index[i] + 1)^-0.1
  #     DF$y[i] <- (DF$x[i] + DF$y[i-1]) * DF$decay[i]
  #   }
  #   
  #   DF
  # }
  # ,
  vector_loop = {
    set.seed(10)
    x <- c(1,sample(c(0:2), replace = T, size = n-1, prob = c(0.8, 0.15, 0.05)))

    y <- vector('numeric', n); y[1] <- 1
    decay <- vector('numeric', n); decay[1] <- 1
    index <- vector('integer', n)
    for(i in 2:n){
      if(x[i] > 0){index[i] <- 0} else {index[i] <- (index[i-1] + 1) }
      decay[i] <- (index[i] + 1)^-0.1
      y[i] <- (x[i] + y[i-1]) * decay[i]
    }
    data.frame(x, y, decay, index)
  }
  ,
  compiled_vector = {
    set.seed(10)
    x <- c(1,sample(c(0:2), replace = T, size = n-1, prob = c(0.8, 0.15, 0.05)))

    base_loop(x)
  }
  ,
  rcpp_func = {
    set.seed(10)
    x <- c(1,sample(c(0:2), replace = T, size = n-1, prob = c(0.8, 0.15, 0.05)))

    decay_func(x)
  }
)

decay_base <- function(x) {


    rle_x <- rle(x > 0)

    index <- sequence(rle_x$lengths)
    index[x != 0] <- 0

    decay <- (index + 1)^(-0.1)

    # initialize y vector and other information
    cum_rle_len <- cumsum(rle_x$lengths)
    y <- vector('numeric', n)
    y[1] <- 1

    # loops through the elements of rle
    for (i in seq_len(length(rle_x$values))[-1]){
      prev_ind <- cum_rle_len[i-1]
      ind_rng <- (prev_ind + 1):(prev_ind + rle_x$lengths[i])

      if (rle_x$values[i]) {
        y[ind_rng] <- y[prev_ind] + cumsum(x[ind_rng])
      } else {
        y[ind_rng] <- cumprod(c(y[prev_ind], decay[ind_rng]))[-1]
      }
    }
    data.frame(x, y, decay, index)

}

base_loop <- function(x) {
  y <- vector('numeric', n); y[1] <- 1
  decay <- vector('numeric', n); decay[1] <- 1
  index <- vector('integer', n)
  for(i in 2:n){
    if(x[i] > 0){index[i] <- 0} else {index[i] <- (index[i-1] + 1) }
    decay[i] <- (index[i] + 1)^-0.1
    y[i] <- (x[i] + y[i-1]) * decay[i]
  }
  data.frame(x, y, decay, index)
}

【讨论】：