将具有大 data.table 的向量合并以执行计算的内存高效方法 (R)

Memory efficient ways of merging a vector with a large data.table to perform calculations (R)

我有一个数据集,其中包含多个模型预测的基于年份的数据,格式为 data.table。

library(data.table)
nYears = 20 # real data: 110
nMod   = 3  # real data: ~ 100
nGrp   = 45

dataset <- data.table(
  group_code = rep(seq(1:nGrp   ),    times= 3*nYears ),
  Year       = rep(seq(1:nYears ),   each=nGrp   ),
  value      = rnorm(2700      , mean = 10, sd = 2),
  var1       = rep (rnorm(nGrp  , mean = nMod, sd = 1) ,  times= nMod*nYears ),
  var2       = rep (rnorm(nGrp  , mean = 1.5, sd = 0.5) , times= nMod*nYears ),
  model   = as.character(rep(seq( from = 1, to = nMod ) ,  each=nGrp  *nYears ))
)
setkey(dataset, Year, model)

我需要基于一个名为 x、长度为 1001 的向量,从这个 数据集 执行一组计算,包含在 seq(-2, 8, by=0.01) 上。 为此,我创建了一个新的 data.table (dt),其中包含重复版本的 dataset 以合并向量 x ,因此:

dt  <- dataset[, lapply(.SD, function(x) rep(x, 1001))]
dt[, x :=  rep(round(seq(-2, 8, by=0.01), 2), each= nYears*nGrp*nMod) ]  

由于我的原始 数据集 包含数百个模型,此操作的内存效率不高。

我需要的最重要的操作包括生成 x 的正态分布,均值 = var1 和 sd= var2,group_code,年份和模型。例如:

 # key computation
 dt [, norm_dist := dnorm (x, var1, var2) , by= .(group_code,  Year, model )]
   

最后这个操作在我的桌面上相当快。但是,我还有其他操作需要执行,这些操作需要对 data.table 进行子集化并且非常消耗 RAM。一个例子:

dt[ x %between% c( 2, 5.99), dt2 := rep_len( rev(dt [x %between% c(-2, 1.99)]$value), length.out=.N) , by= .(Year, model) ] 

弹出以下错误:

Error: cannot allocate vector of size 1.3 Gb

我认为此特定步骤中的问题与子集和 rev() 函数有关。

然而,我使用的基于来自 data.table dt 的向量“x”执行一组计算的方法似乎并不合适,因为我将数据集与计算所需的向量合并 ("x").

我希望有人能教我如何有效地改进我的代码,因为我在原始数据集中有相当多的模型,大大增加了它的大小。

谢谢!

我觉得这部分代码应该更清晰

dt[ x %between% c( 2, 5.99), dt2 := rep_len( rev(dt [x %between% c(-2, 1.99)]$value), length.out=.N) , by= .(Year, model) ]

因为它对我来说有点像黑盒子。特别是因为这个双重子集是你的问题产生的地方。

这些代码位,x %between% c( 2, 5.99)dt[x %between% c(-2, 1.99)],在您的所有情况下应该始终位于相同的位置。您应该在代码中考虑这一点以提高效率。

尝试这样的事情让事情更清楚一点:

by_YM <- split(dt, by=c("Year", "model"))
ind1  <- which(by_YM[[1]][["x"]] %between% c( 2, 5.99))
ind2  <- which(by_YM[[1]][["x"]] %between% c(-2, 1.99))

for(i in 1:length(by_YM)){
 
  dt_i <- by_YM[[i]]
  #val1 <- rep_len(rev(dt_i$value[ind2]), length.out=length(ind1)) #val1 is equal to val, no need for rep_len
  val  <- rev(dt_i$value[ind2])
  
  by_YM[[i]] <- dt_i[ind1, dt2 := val] 
  
}

但是我们的 dt2 列不相等,但由于我不确定最终结果应该如何,我无法调试它。

dt2_a <- dt[Year == 20 & model == 3, dt2]
dt2_b <- by_YM[["20.3"]][, dt2]

test  <- cbind(dt2_a, dt2_b)

第二个代码也快得多。

library(microbenchmark)

microbenchmark( "new_code" = {
  by_YM <- split(dt, by=c("Year", "model"))

ind1  <- which(by_YM[[1]][["x"]] %between% c( 2, 5.99))
ind2  <- which(by_YM[[1]][["x"]] %between% c(-2, 1.99))

for(i in 1:length(by_YM)){
  
  dt_i <- by_YM[[i]]
  val1 <- rep_len(rev(dt_i$value[ind2]), length.out=length(ind1)) #val1 is equal to val, no need for rep_len
  val  <- rev(dt_i$value[ind2])
  
  by_YM[[i]] <- dt_i[ind1, dt2 := val] 
  
}}, "old_code" = dt[ x %between% c( 2, 5.99), 
                   dt2 := rep_len( rev(dt [x %between% c(-2, 1.99)]$value), length.out=.N) , by= .(Year, model) ],  
times = 5)

Unit: milliseconds
     expr      min        lq      mean    median        uq       max neval cld
 new_code  155.426  156.4916  200.6587  185.0347  188.9436  317.3977     5  a 
 old_code 1290.909 1299.8570 1398.6866 1370.4526 1471.0569 1561.1574     5   b

试试吧,祝你好运