作者热门文章
- html - 出于某种原因,IE8 对我的 Sass 文件中继承的 html5 CSS 不友好?
- JMeter 在响应断言中使用 span 标签的问题
- html - 在 :hover and :active? 上具有不同效果的 CSS 动画
- html - 相对于居中的 html 内容固定的 CSS 重复背景?
顺序累积计算
我需要进行时间序列计算,其中每一行计算的值取决于前一行计算的结果。我希望使用 data.table
的便利性.实际问题是一个水文模型——累积水平衡计算,在每个时间步增加降雨量并减去作为当前水量函数的径流和蒸发量。数据集包括不同的盆地和情景(组)。在这里,我将使用一个更简单的问题来说明问题。
计算的简化示例如下所示,对于每个时间步(行)i
:
v[i] <- a[i] + b[i] * v[i-1]
a
和
b
是参数值的向量,
v
是结果向量。对于第一行 (
i == 1
) 的初始值
v
被视为
v0 = 0
.
shift()
在
data.table
.一个最小的例子,包括想要的结果
v.ans
, 是
library(data.table) # version 1.9.7
DT <- data.table(a = 1:4,
b = 0.1,
v.ans = c(1, 2.1, 3.21, 4.321) )
DT
# a b v.ans
# 1: 1 0.1 1.000
# 2: 2 0.1 2.100
# 3: 3 0.1 3.210
# 4: 4 0.1 4.321
DT[, v := NA] # initialize v
DT[, v := a + b * ifelse(is.na(shift(v)), 0, shift(v))][]
# a b v.ans v
# 1: 1 0.1 1.000 1
# 2: 2 0.1 2.100 2
# 3: 3 0.1 3.210 3
# 4: 4 0.1 4.321 4
这不起作用,因为
shift(v)
提供原始列
v
的副本,移动了 1 行。它不受分配给
v
的影响.
vcalc <- function(a, b, v0 = 0) {
v <- rep(NA, length(a)) # initialize v
for (i in 1:length(a)) {
v[i] <- a[i] + b[i] * ifelse(i==1, v0, v[i-1])
}
return(v)
}
此累积函数适用于 data.table:
DT[, v := vcalc(a, b, 0)][]
# a b v.ans v
# 1: 1 0.1 1.000 1.000
# 2: 2 0.1 2.100 2.100
# 3: 3 0.1 3.210 3.210
# 4: 4 0.1 4.321 4.321
identical(DT$v, DT$v.ans)
# [1] TRUE
我的问题
data.table
方式,而不必使用 for 循环和/或函数定义?使用
set()
也许?
ifelse()
来自
for
环形:
vcalc2 <- function(a, b, v0 = 0) {
v <- rep(NA, length(a))
for (i in 1:length(a)) {
v0 <- v[i] <- a[i] + b[i] * v0
}
return(v)
}
vcalc2()
比
vcalc()
快 60% .
最佳答案
它可能不是您正在寻找的 100%,因为它不使用“data.table-way”并且仍然使用 for 循环。但是,这种方法应该更快(我假设您想使用 data.table 和 data.table-way 来加速您的代码)。我利用 Rcpp 编写了一个名为 HydroFun
的简短函数。 ,可以像任何其他函数一样在 R 中使用(您只需要先获取该函数)。我的直觉告诉我 data.table 方式(如果存在)非常复杂,因为您无法计算封闭形式的解决方案(但我在这一点上可能是错误的......)。
我的方法如下所示:
Rcpp 函数如下所示(在文件中:hydrofun.cpp
):
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
NumericVector HydroFun(NumericVector a, NumericVector b, double v0 = 0.0) {
// get the size of the vectors
int vecSize = a.length();
// initialize a numeric vector "v" (for the result)
NumericVector v(vecSize);
// compute v_0
v[0] = a[0] + b[0] * v0;
// loop through the vector and compute the new value
for (int i = 1; i < vecSize; ++i) {
v[i] = a[i] + b[i] * v[i - 1];
}
return v;
}
Rcpp::sourceCpp("hydrofun.cpp")
library(data.table)
DT <- data.table(a = 1:4,
b = 0.1,
v.ans = c(1, 2.1, 3.21, 4.321))
DT[, v_ans2 := HydroFun(a, b, 0)]
DT
# a b v.ans v_ans2
# 1: 1 0.1 1.000 1.000
# 2: 2 0.1 2.100 2.100
# 3: 3 0.1 3.210 3.210
# 4: 4 0.1 4.321 4.321
library(microbenchmark)
n <- 10000
dt <- data.table(a = 1:n,
b = rnorm(n))
microbenchmark(dt[, v1 := vcalc(a, b, 0)],
dt[, v2 := HydroFun(a, b, 0)])
# Unit: microseconds
# expr min lq mean median uq max neval
# dt[, `:=`(v1, vcalc(a, b, 0))] 28369.672 30203.398 31883.9872 31651.566 32646.8780 68727.433 100
# dt[, `:=`(v2, HydroFun(a, b, 0))] 381.307 421.697 512.2957 512.717 560.8585 1496.297 100
identical(dt$v1, dt$v2)
# [1] TRUE
关于r - 如何在 data.table 中编写累积计算,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/40412516/
我是一名优秀的程序员,十分优秀!