确定动态 window 宽度：满足条件的值的有效滚动计数

Question

我有一个包含两列 a 和 b 的 data.frame，其中 a 已排序。我想获得 b 的滚动平均值，其中 window 是 a - 5 到 a 的范围（即从 a 的当前值到 a - 5 是).

使用 data.table::frollmean() 执行具有不同 window 宽度的滚动平均值是微不足道的（adaptive = TRUE；“每个单独的观察值都有自己相应的滚动 window 宽度”），所以唯一的问题是计算那些 window 宽度。

那么，给定以下 data.frame，我如何确定每个均值的 window 大小？

set.seed(42)
x <- data.frame(
    a = sort(runif(10, 0, 10)),
    b = 1:10
)
x
#>           a  b
#> 1  1.346666  1
#> 2  2.861395  2
#> 3  5.190959  3
#> 4  6.417455  4
#> 5  6.569923  5
#> 6  7.050648  6
#> 7  7.365883  7
#> 8  8.304476  8
#> 9  9.148060  9
#> 10 9.370754 10

^{由 reprex package (v0.3.0)}

于 2020-07-03 创建

如果我将 window 大小作为新列 n，我希望结果是

#>           a  b n
#> 1  1.346666  1 1
#> 2  2.861395  2 2
#> 3  5.190959  3 3
#> 4  6.417455  4 3
#> 5  6.569923  5 4
#> 6  7.050648  6 5
#> 7  7.365883  7 6
#> 8  8.304476  8 6
#> 9  9.148060  9 7
#> 10 9.370754 10 8

所以，比如a[2] = 2.86和2.86 - 5之间有两个值（包括它本身），a[8] = 8.30和8.30 - 5之间有六个值。

我已经使用 outer 做到了：

suppressPackageStartupMessages({
    library(magrittr)
    library(data.table)
})

f <- function(x, y) {
    return(y %between% list(x - 5, x))
}

outer(x$a, x$a, f) %>% rowSums()
#>  [1] 1 2 3 3 4 5 6 6 7 8

然而，我的实际案例有 5000 行，而且这种方法非常慢（大约需要 10 秒）。我看到的一个问题是它将 a 的每个值与 a 的每个其他值进行比较，因此必须执行大约 25,000,000 次比较。但是，我知道 a 是排序的，所以如果我们在比较中找到一段 TRUE 结果，然后是 FALSE，我们知道当前值 [=14] 的所有后续结果=] 也将是 FALSE（这意味着我们在允许范围内，然后超过了最高允许值 a，因此其他所有内容也将被拒绝）。

那么，有更好、更快的方法吗？

Answer 1

因为您似乎无论如何都会加载 data.table（对于 frollmean），您可以将 data.frame 强制为 data.table，并通过引用添加新列.

findInterval用于查找每个减去的值在原始值中的索引。然后从通过 .I 或 seq_along 获得的原始索引中减去该索引，以获得 window 大小。

setDT(x)
x[ , n := .I - findInterval(a - 5, a)]

# x
#            a  b n
#  1: 1.346666  1 1
#  2: 2.861395  2 2
#  3: 5.190959  3 3
#  4: 6.417455  4 3
#  5: 6.569923  5 4
#  6: 7.050648  6 5
#  7: 7.365883  7 6
#  8: 8.304476  8 6
#  9: 9.148060  9 7
# 10: 9.370754 10 8

类似于base：

x$n = seq_along(x$a) - findInterval(x$a - 5, x$a)

Answer 2

这是另一种方法，在非 equi 自连接中聚合：

library(data.table)
setDT(x)[, low := a - 5][
  , n := x[x, on = .(a >= low , a <= a), by = .EACHI, .N]$N][
      , low := NULL][]

           a  b n
 1: 1.346666  1 1
 2: 2.861395  2 2
 3: 5.190959  3 3
 4: 6.417455  4 3
 5: 6.569923  5 4
 6: 7.050648  6 5
 7: 7.365883  7 6
 8: 8.304476  8 6
 9: 9.148060  9 7
10: 9.370754 10 8

但是 OP 的目标是 计算具有变量 window 大小 的滚动平均值。

所以，既然我们可以一次性搞定，为什么要停下来调用 frollmean()？:

library(data.table)
setDT(x)[, low := a - 5][
  , roll.mean := x[x, on = .(a >= low , a <= a), by = .EACHI, mean(b)]$V1][
    , low := NULL][]

           a  b roll.mean
 1: 1.346666  1       1.0
 2: 2.861395  2       1.5
 3: 5.190959  3       2.0
 4: 6.417455  4       3.0
 5: 6.569923  5       3.5
 6: 7.050648  6       4.0
 7: 7.365883  7       4.5
 8: 8.304476  8       5.5
 9: 9.148060  9       6.0
10: 9.370754 10       6.5

基准

由于 OP 关心他的生产用例的性能，这里是一个基准，它改变了行数以及 window:

的大小

library(bench)
library(ggplot2)

bm <- press(
  n = 10^(c(2, 3, 4)),
  window_size = c(5, 15, 50),
  {
    set.seed(42)
    x0 <- data.table(
      a = sort(runif(n, 0, n)),
      b = seq(n)
    )
    mark(
      findInterval = {
        x <- copy(x0)
        x[, roll.mean := frollmean(b, .I - findInterval(a - window_size, a), adaptive = TRUE)]
      },
      non_equi_join = {
        x <- copy(x0)
        x[, low := a - window_size][
          , roll.mean := x[x, on = .(a >= low , a <= a), by = .EACHI, mean(b)]$V1][
            , low := NULL]
      }
    )
  }
)

autoplot(bm)

显然，

与自适应 frollmean() 的组合总是比 non-equi join 方法
window大小似乎对性能没有影响。

确定动态 window 宽度：满足条件的值的有效滚动计数

Determine dynamic window width: efficient rolling count of values which satisfy a condition

r

data.table

rolling-computation

基准