r - 在 for 循环中对 data.table 进行子集化较慢且资源匮乏

Question

使用 data.table R 包时，我注意到在运行一个简单的 for 循环时处理器使用率非常高，该循环将使用来自另一个 data.table 的值对数据集进行子集化。当我说高使用率时，我的意思是在循环运行的整个时间内 100% 的所有可用线程。

有趣的部分是，对相同的进程使用 data.frame 对象对相同的输出花费的时间少 10 倍。并且只有一个核心达到 100%。

这是我希望可重现的示例:

chr = c(rep(1, 1000), rep(2, 1000), rep(3, 1000), rep(3,1000))
start = rep(seq(from =1, to = 100000, by=100), 4)
end = start + 100

df1 <- data.frame(chr=chr, start=start, end=end)
df2 <- rbind(df1,df1,df1,df1,df1)
dt1 <- data.table::data.table(df1)
dt2 <- data.table::data.table(df2)

test1 <- list()
test2 <- list()

#loop subsetting a data.frame
system.time(
for (i in 1:nrow(df2)) {
  no.dim <- dim(df1[df1$chr == df2[i, 'chr'] & df1$start >= df2[i, 'start'] & df1$end <= df2[i, 'end'], ])[1]
  test1[i] <- no.dim
})

# loop subsetting a data.table using data.table syntax
system.time(
for (i in 1:nrow(dt2)) {
  no.dim <- dim(dt1[chr == dt2[i, chr] & start >= dt2[i, start] & end <= dt2[i, end], ])[1]
  test2[i] <- no.dim
})

# is the output the same
identical(test1, test2)

这是输出:

> #loop subsetting a data.frame
> system.time(
+ for (i in 1:nrow(df2)) {
+   no.dim <- dim(df1[df1$chr == df2[i, 'chr'] & df1$start >= df2[i, 'start'] & df1$end <= df2[i, 'end'], ])[1]
+   test1[i] <- no.dim
+ })
   user  system elapsed 
  2.607   0.004   2.612 
> 
> # loop subsetting a data.table using data.table syntax
> system.time(
+ for (i in 1:nrow(dt2)) {
+   no.dim <- dim(dt1[chr == dt2[i, chr] & start >= dt2[i, start] & end <= dt2[i, end], ])[1]
+   test2[i] <- no.dim
+ })
   user  system elapsed 
192.632   0.152  24.398 
> 
> # is the output the same
> identical(test1, test2)
[1] TRUE

现在，我知道可能有多种更好、更有效的方法来执行相同的任务，而且我可能不是按照 data.table 的方式来做的。但是假设出于某种原因，您有一个使用“data.frame”对象的脚本，并且您想快速重写该内容以改用 data.table。上面采取的方法似乎完全有道理。

任何人都可以重现关于减速和高处理器使用率的相同情况吗？它是否可以通过保持或多或少相同的子集化过程来修复，还是必须完全重写才能在 data.table 上有效使用？

PS:刚在windows机器上测试过，线程使用正常(一个线程100%运行)，但还是比较慢。在类似于我的系统上测试它给出了与上面相同的结果。

R version 3.5.1 (2018-07-02)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 18.10

Matrix products: default
BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.8.0
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.8.0

locale:
 [1] LC_CTYPE=C           LC_NUMERIC=C         LC_TIME=C            LC_COLLATE=C        
 [5] LC_MONETARY=C        LC_MESSAGES=C        LC_PAPER=et_EE.UTF-8 LC_NAME=C           
 [9] LC_ADDRESS=C         LC_TELEPHONE=C       LC_MEASUREMENT=C     LC_IDENTIFICATION=C 

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] data.table_1.12.0

loaded via a namespace (and not attached):
 [1] compiler_3.5.1   assertthat_0.2.0 cli_1.0.1        tools_3.5.1      pillar_1.3.1    
 [6] rstudioapi_0.9.0 tibble_2.0.0     crayon_1.3.4     utf8_1.1.4       fansi_0.4.0     
[11] pkgconfig_2.0.2  rlang_0.3.1   

编辑:

感谢大家的评论。看来减速问题与@Hugh 详述的 [.data.table 的开销有关。正如@denis 所指出的，在 efficient subsetting of data.table with greater-than, less-than using indices 中提到了同样的问题。

@Frank 提出的修复虽然确实有效并产生类似的输出，但通过完全删除循环并在原始数据集中添加可能不需要的列来改变过程的行为。

编辑.1:

在我第一次编辑后，@Frank 添加了另一种方法，该方法包括使用 data.table 语法计算列表列。虽然它非常整洁，但我必须承认我需要一段时间来弄清楚发生了什么。我认为它只是在子集 data.table 的开始和结束列上计算 lm()，所以我尝试使用 for 循环和 data.frames 重现结果。时间:

> system.time({res <- dt1[dt2, on=.(chr, start >= start, end <= end), .(n = .N, my_lm = list(lm(x.start ~ x.end))), by=.EACHI][, .(n, my_lm)]; res <- as.list(res$my_lm)})
   user  system elapsed 
 11.538   0.003  11.336 
> 
> test_new <- list()
> system.time(
+   for (i in 1:20000) {
+     df_new <- df1[df1$chr == df2$chr[i] & df1$start >= df2$start[i] & df1$end <= df2$end[i],]
+     test_new[[i]] <- lm(df_new$start ~ df_new$end)
+   })
   user  system elapsed 
 12.377   0.048  12.425 
> 

只要你有 lm() 这样的瓶颈函数，你最好使用基本的 for 循环(为了控制和可读性)，但使用 data.frames。

Answer 1

用户时间和耗时之间的差异是一个线索，表明在幕后进行了一些并行化:

library(data.table)
chr = c(rep(1, 1000), rep(2, 1000), rep(3, 1000), rep(3,1000))
start = rep(seq(from =1, to = 100000, by=100), 4)
end = start + 100

df1 <- data.frame(chr=chr, start=start, end=end)
df2 <- rbind(df1,df1,df1,df1,df1)
dt1 <- data.table::data.table(df1)
dt2 <- data.table::data.table(df2)

print(dim(dt1))
#> [1] 4000    3
print(dim(dt2))
#> [1] 20000     3


test1 <- list()
test2 <- list()

bench::system_time({
  for (i in 1:nrow(df2)) {
    no.dim <- dim(df1[df1$chr == df2[i, 'chr'] &
                        df1$start >= df2[i, 'start'] &
                        df1$end <= df2[i, 'end'], ])[1]
    test1[i] <- no.dim
  }
})
#> process    real 
#>  3.547s  3.549s

print(getDTthreads())
#> [1] 12

bench::system_time({
  for (i in 1:nrow(dt2)) {
    no.dim <- dim(dt1[chr == dt2[i, chr] & start >= dt2[i, start] & end <= dt2[i, end], ])[1]
    test2[i] <- no.dim
  }
})
#> process    real 
#> 83.984s 52.266s

setDTthreads(1L)
bench::system_time({
  for (i in 1:nrow(dt2)) {
    no.dim <- dim(dt1[chr == dt2[i, chr] & start >= dt2[i, start] & end <= dt2[i, end], ])[1]
    test2[i] <- no.dim
  }
})
#> process    real 
#> 30.922s 30.920s

由 reprex package (v0.2.1) 于 2019-01-30 创建

但同样重要的是您调用 [ 20,000 次。考虑这个最小的用途来证明单行表的 [.data.table 的开销在运行时占主导地位:

library(data.table)
chr = c(rep(1, 1000), rep(2, 1000), rep(3, 1000), rep(3,1000))
start = rep(seq(from =1, to = 100000, by=100), 4)
end = start + 100

df1 <- data.frame(chr=chr, start=start, end=end)
df2 <- rbind(df1,df1,df1,df1,df1)
dt1 <- data.table::data.table(df1)
dt2 <- data.table::data.table(df2)

bench::system_time({
  o <- integer(nrow(df2))
  for (i in 1:nrow(df2)) {
    o[i] <- df2[i, ][[2]]
  }
})
#>   process      real 
#> 875.000ms 879.398ms

bench::system_time({
  o <- integer(nrow(dt2))
  for (i in 1:nrow(dt2)) {
    o[i] <- dt2[i, ][[2]]
  }
})
#> process    real 
#> 26.219s 13.525s

由 reprex package (v0.2.1) 于 2019-01-30 创建