提高 R 的效率（矢量化？）

Question

我在 R 中有一个脚本需要 8 分钟才能运行，它基本上比较了多年期间 800 条记录的日期范围。这太长了。我是 R 新手，很确定它与我的嵌入式循环有关。此外，当我尝试将我的数据转换为玩具问题时，它似乎不起作用。我一直在处理从 excel 中读取的数组类型。

# data vectors
ID <- c("1e", "1f", "1g")
StartDate <- c(1, 2, 4)
EndDate <- c(3, 4, 5)
Type <- c("A", "B", "B")
Qty <- c(.5, 2.5, 1)

# table rows and headers
Days <- c(1, 2, 3, 4, 5)
setOfTypes <- c("A", "B")

# get subset of active IDs for each day in table
ActiveID <- data.frame()
for(d in 1:length(Days)){
  check <- StartDate<=Days[d] & EndDate>=Days[d]
  subsetID <- subset(ID, check)
  strSubsetID <- c()
  for(i in 1:length(subsetID)){
    strSubsetID <- paste(ID, subsetID[i], sep=",")
}
ActiveID[d,1] <- strSubsetID
}

# calculate quantity counts by day and type
Count <- matrix(,length(Days),length(setOfTypes))
for(d in 1:length(Days)){
  for(t in 1:length(setOfTypes))
    check <- Type == setOfTypes[t] & sapply(ID, grepl, x=ActiveID[d,1])
    tempCount <- subset(Types, check)
    Count[t,d] <- sum(tempCount)
  }
}

结果应该是一个表格（天 x 类型），其中每个元素由给定日期和类型的活动 ID 的数量总和组成。

我希望将此代码矢量化，以便在应用于更大的数据集时运行得更快！！请帮忙，谢谢。

Answer 1

您的代码没有按原样运行，因此我无法确切知道您在寻找什么。您的描述表明您希望 Days 中的每个 StartDate 和 EndDate 之间的 Qty 总和，按 Type 分组。这将产生这样一个矩阵：

df <- data.frame(ID,StartDate,EndDate,Type,Qty,stringsAsFactors=FALSE)
Days <- min(StartDate):max(EndDate)

is.between <- function(x,df) with(df,x>=StartDate & x<=EndDate)
get.sums   <- function(df) sapply(Days,function(d,df) sum(df[is.between(d,df),"Qty"]),df)
do.call(rbind,lapply(split(df,df$Type), get.sums))
#   [,1] [,2] [,3] [,4] [,5]
# A  0.5  0.5  0.5  0.0    0
# B  0.0  2.5  2.5  3.5    1

这是一种可能更快的 data.table 方法。注意is.between(...) 和get.sums(...) 的不同定义。

DT <- data.table(df,key="Type")
is.between <- function(x,a,b) x>=a & x <= b
get.sums   <- function(day) DT[,list(day,Qty=sum(Qty[is.between(day,StartDate,EndDate)])),by=Type]
long       <- rbindlist(lapply(Days,get.sums))
result     <- dcast.data.table(long,Type~day,value.var="Qty")
result
#    Type   1   2   3   4 5
# 1:    A 0.5 0.5 0.5 0.0 0
# 2:    B 0.0 2.5 2.5 3.5 1

---

这里有一些基准测试，其中有希望更具代表性的示例数据集（800 行，500 个开始日期，总日期范围 >900 天），并且还测试了 @Arun 的答案。

# more representative example
set.seed(1)  # for reproducibility
StartDate <- sample(1:500,800,replace=TRUE)
EndDate   <- StartDate + rpois(800,400)
Type      <- sample(LETTERS[1:20],800,replace=TRUE)
Qty       <- rnorm(800,10,2)
Days      <- min(StartDate):max(EndDate)
df        <- data.frame(StartDate,EndDate,Type,Qty, stringsAsFactors=FALSE)

数据框方法与两种数据表方法的比较。

library(data.table)
library(reshape2)
DT <- data.table(df,key="Type")
f.df <- function() {
  is.between <- function(x,df) with(df,x>=StartDate & x<=EndDate)
  get.sums   <- function(df) sapply(Days,function(d,df) sum(df[is.between(d,df),"Qty"]),df)
  do.call(rbind,lapply(split(df,df$Type), get.sums))
}
f.dt1 <- function() {
  is.between <- function(x,a,b) x>=a & x <= b
  get.sums   <- function(day) DT[,list(day,Qty=sum(Qty[is.between(day,StartDate,EndDate)])),by=Type]
  long       <- rbindlist(lapply(Days,get.sums))
  dcast.data.table(long,Type~day,value.var="Qty")
}
f.dt2 <- function() {
  lookup <- data.table(StartDate=Days, EndDate=Days)
  setkey(lookup)
  j_olaps <- foverlaps(DT, lookup, by.x=c("StartDate", "EndDate"), type="any")
  dcast.data.table(j_olaps, Type ~ StartDate, value.var="Qty", fun.agg=sum, na.rm=TRUE)
}
identical(f.dt1(),f.dt2())   # same result? YES!
# [1] TRUE
library(microbenchmark)
microbenchmark(f.df(),f.dt1(),f.dt2(),times=10)
# Unit: milliseconds
#     expr        min         lq    median        uq       max neval
#   f.df() 1199.76370 1212.03787 1222.6558 1243.8743 1275.5526    10
#  f.dt1() 1634.92675 1664.98885 1689.7812 1714.2662 1798.9121    10
#  f.dt2()   91.53245   95.19545  129.2789  158.0789  208.1818    10

所以@Arun 的方法比 df 方法快约 10 倍，比上面的 dt 方法快约 17 倍。

Answer 2

查看@jihoward 的代码，这似乎是重叠连接 的一个案例，它最近在data.table 的v1.9.4 中实现。该函数称为foverlaps()。下面是我们如何使用它：

首先，我们创建一个查找表，其中包含我们希望重叠连接的日期范围。这是使用@jihoward 代码中的变量Days 构造的。在您的情况下，开始日期和结束日期是相同的。

require(data.table) ## 1.9.4
lookup <- data.table(StartDate=Days, EndDate=Days)
setkey(lookup)

然后我们使用foverlaps() 计算重叠连接。这里的重叠类型指定为any。查看?foverlaps 以了解这意味着什么，以及可以执行的其他类型的重叠。

j_olaps = foverlaps(DT, lookup, by.x=c("StartDate", "EndDate"), type="any")

现在我们已经有了重叠，我们可以简单地将其转换如下：

dcast.data.table(j_olaps, Type ~ StartDate, value.var="Qty", fun.agg=sum, na.rm=TRUE)

#    Type   1   2   3   4 5
# 1:    A 0.5 0.5 0.5 0.0 0
# 2:    B 0.0 2.5 2.5 3.5 1

我相信这应该比必须对Days 中的每个元素进行基于矢量扫描的子集要快得多。很高兴知道您获得了多少加速（如果有的话）！

HTH