不是答案,但可能有助于解决问题。似乎最坏情况下的性能是对许多短组求和,这似乎与向量的大小成线性关系
> n = 100000; x = runif(n); f = sample(n/2, n, TRUE)
> system.time(rowsum(x, f))
user system elapsed
0.228 0.000 0.229
> n = 1000000; x = runif(n); f = sample(n/2, n, TRUE)
> system.time(rowsum(x, f))
user system elapsed
1.468 0.040 1.514
> n = 10000000; x = runif(n); f = sample(n/2, n, TRUE)
> system.time(rowsum(x, f))
user system elapsed
17.369 0.748 18.166
似乎有两个捷径可用,避免重新订购
> n = 10000000; x = runif(n); f = sample(n/2, n, TRUE)
> system.time(rowsum(x, f, reorder=FALSE))
user system elapsed
16.501 0.476 17.025
并避免对性格的内在胁迫
> n = 10000000; x = runif(n); f = as.character(sample(n/2, n, TRUE));
> system.time(rowsum(x, f, reorder=FALSE))
user system elapsed
8.652 0.268 8.949
然后似乎涉及的基本操作 - 找出分组因子的唯一值(预先分配结果向量)并求和
> n = 10000000; x = runif(n); f = sample(n/2, n, TRUE)
> system.time({ t = tabulate(f); sum(x) })
user system elapsed
0.640 0.000 0.643
所以是的,似乎有相当大的空间可以实现更快的单一用途。这似乎是一个自然的 for data.table
,并且在 C 中实现起来并不难。这是一个混合解决方案,使用 R 进行制表,使用“经典”C 接口进行求和
library(inline)
rowsum1.1 <- function(x, f) {
t <- tabulate(f)
crowsum1(x, f, t)
}
crowsum1 = cfunction(c(x_in="numeric", f_in="integer", t_in = "integer"), "
SEXP res_out;
double *x = REAL(x_in), *res;
int len = Rf_length(x_in), *f = INTEGER(f_in);
res_out = PROTECT(Rf_allocVector(REALSXP, Rf_length(t_in)));
res = REAL(res_out);
memset(res, 0, Rf_length(t_in) * sizeof(double));
for (int i = 0; i < len; ++i)
res[f[i] - 1] += x[i];
UNPROTECT(1);
return res_out;
")
和
> system.time(r1.1 <- rowsum1.1(x, f))
user system elapsed
1.276 0.092 1.373
要实际返回与 相同的结果,rowsum
需要将其成形为具有适当暗名称的矩阵
rowsum1 <- function(x, f) {
t <- tabulate(f)
r <- crowsum1(x, f, t)
keep <- which(t != 0)
matrix(r[keep], ncol=1, dimnames=list(keep, NULL))
}
> system.time(r1 <- rowsum1(x, f))
user system elapsed
9.312 0.300 9.641
所以对于所有这些工作,我们的速度只有 2 倍(而且更不通用——x 必须是数字,f 必须是整数;没有 NA 值)。是的,存在效率低下,例如,分配没有计数的空间级别(尽管这避免了对名称的字符向量进行昂贵的强制)。