从 data.frame 中的每个记录时间戳计算一秒内出现多少 data.frame 记录时间戳的优化 window

Optimizaion of calculating how many data.frame record timestamps appear in a one second window from every record timestamp in that data.frame

我有一个包含多个列的 R data.frame,其中一列包含 POSIXct 时间戳记记录。我想在 data.frame 中添加一列,对于每一行,该列包含时间戳介于该行的时间戳和未来一秒之间的记录数。

下面的代码实现了这一点,但它对我正在处理的数据(通常是 60K+ 条记录)来说真的很慢。我想知道是否有办法加快速度。

# Create a data frame with POSIXct values spread over a few minutes.
# The actual number of records can be over 60,000.
set.seed(1234)
times <- as.POSIXct("2015-02-18 11:39:17.104206 AEDT") + 
    runif(n = 10000, min = -5*60, max = 5*60)
times <- sort(times) # my source data comes to me sorted
times <- data.frame(datetime = times)

# For each event (timestamp), calculate how many events (timestamps) appear in
# a one second window following that event.
system.time(
for (i in 1:length(times$datetime)) {
        times$eventCount[i] <- sum(
                times$datetime >= times$datetime[i] & 
                times$datetime < times$datetime[i] + 1)
}
)

我系统上的结果是:

user  system elapsed
8.10    0.00    8.21

有趣的是,处理时间与记录数不成线性关系。对于 20K 条记录,用户时间为 24.74 秒。

查看类似的问题(如 this 一个和其中引用的问题)会建议使用 data.table 应该大大加快速度,但我无法弥合代码之间的差距在这些答案中(查看给定记录两侧的固定数量的记录)到我需要的(查看给定记录一侧的未知数量的记录)。

rcpp 看起来是最好的方法,但我根本不懂任何 c++。

感谢任何帮助!

当 hte 数量增加时,它对我来说增加了一倍以上。我虽然也许我可以通过避免使用无效的“$”访问数据框来获得更好的性能,但我确实找到了一种节省时间的方法。您不应该测试少于当前时间的次数,因为您知道 (i-1) 的答案已经给出了数据的排序性质。只需记录 1 秒内前面的项目数。 (我的处理器和你的差不多,所以这个结果实际上比第一个快了大约 25%:

system.time( {dt <- times$datetime
 for (i in 1:length(dt)) {
         eventCount[i] <- sum(
                 tail(dt, 10000-i) < dt[i] + 1)
 }}
 )
   user  system elapsed 
  5.410   0.716   6.042 

尝试

library(dplyr)
dt <- times$datetime
system.time({
times <- times %>% 
  rowwise() %>% 
  mutate(eventCount = sum(between(dt, datetime, datetime+1))) 
})

在 AWS 免费套餐上,

user  system elapsed 
3.309   0.048   3.358 

[编辑]

显然,dplyr 中的 between 相当慢。在 Rcpp 中实现这一步很容易,结果非常有成果。 betweenCpp(v,x,y) 的行为类似于 sum(between(v, x, y)),即计算 'v' 中位于 xy 之间的元素。

Rcpp::cppFunction('int betweenCpp(NumericVector v, double x, double y) {
                  NumericVector::iterator low1, low2;
                  low1=std::lower_bound (v.begin(), v.end(), x);
                  low2=std::lower_bound (v.begin(), v.end(), y);
                  return (low2 - low1);
                  }')

现在,我们可以使用 Rcpp 的强大功能,在我的普通笔记本电脑上 n=100.000 花费了 0.28 秒。

dt <- times$datetime
system.time({
  times <- times %>% 
    rowwise() %>% 
    mutate(eventC = betweenCpp(dt, datetime, datetime+1)) 
})

[附录]

如果你想要额外的速度,你可以做完整的 Rcpp 实现。

#include <Rcpp.h>
#include <algorithm>
using namespace Rcpp;

int betweenCpp(NumericVector v, double x, double y) {
  NumericVector::iterator low1, low2;
  low1=std::lower_bound (v.begin(), v.end(), x);
  low2=std::lower_bound (v.begin(), v.end(), y);
  return (low2 - low1);
}
// [[Rcpp::export]]
NumericVector EventCountCpp(NumericVector x) {
  int n=x.size();
  NumericVector count(n);
  for (int i = 0; i < n; i++) {
    count[i]=betweenCpp(x, x[i], x[i]+1);
  }
  return(count);
}

在您的工作目录中将其保存为 count.cpp,然后 n=100.000.

花费了 0.01 秒
sourceCpp("count.cpp")
system.time(times$EventCount <- EventCountCpp(times$datetime))

基于@Kashaa 的 Rcpp 解决方案的更简单的逻辑。

数据

require(dplyr)
require(data.table)
set.seed(1234L)
dt = data.table(datetime=as.POSIXct("2015-02-18 11:39:17.104206 AEDT") + 
    runif(n = 100000, min = -5*60, max = 5*60), key="datetime")
df = as.data.frame(dt)

data.table解决方案

setNumericRounding(0L)
betweendt <- function(x, col, eps) {
    idx1 = dt[.(col), mult="first", roll=-Inf, which=TRUE]
    idx2 = dt[.(col+1-eps-unclass(col)*eps), 
                mult="last", roll=Inf, which=TRUE]
    idx2-idx1+1L
}
system.time({
dt[, eventC := betweendt(dt, dt$datetime, .Machine$double.eps)]
})
#    user  system elapsed 
#   0.043   0.001   0.045 

Rcpp 版本(来自@Khashaa)

system.time({
  col = df$datetime
  df <- df %>% 
    rowwise() %>% 
    mutate(eventC = betweenCpp(col, datetime, datetime+1)) 
})
#    user  system elapsed 
#   0.142   0.001   0.142 

identical(df$eventC, dt$eventC)
# [1] TRUE

data.table 这里的解决方案快了约 3 倍。


使用 foverlaps() 参考旧版本的历史记录(这太过分了)。