将数字向量与另一个向量中的值进行比较的最佳方法

Question

假设我有两个值向量：

a <- c(1,3,4,5,6,7,3)
b <- c(3,5,1,3,2)

而且我想将一些函数 ,FUN应用于 的每个输入，以a对抗整个b，最有效的方法是什么。

更具体地说，在这种情况下，对于a'a' 的每个值，对于每个元素，我想知道有多少元素b大于或等于该值。天真的方法是执行以下操作：

sum(a < b)

当然，这不起作用，因为它试图并行迭代每个向量并给我警告：

较长的物体长度不是较短物体长度的倍数

顺便说一句，该命令的输出是3.

但是，在我的情况下，我希望看到的输出是：

0 2 4 4 5 5 2

当然，我意识到我可以使用 for 循环来做到这一点：

out <- c()
for (i in a) {
    for (i in a) { out[length(out) + 1] = sum(b<i)}
}

同样，我可以这样使用sapply：

sapply(a, function(x)sum(b<x))

但是，我正在努力成为一名优秀的 R 程序员并远离 for 循环并且sapply似乎很慢。还有其他选择吗？

对于它的价值，我这样做了几百万次，其中length(b)总是小于1 到 30 length(a)，length(a)范围是 1 到 30。

Answer 1

尝试这个：

findInterval(a - 0.5, sort(b))

从 a) 避免sort和 b) 避免开销findInterval和order使用更简单的.Internal包装器来提高速度：

order2 = function(x) .Internal(order(T, F, x))

findInterval2 = function(x, vec, rightmost.closed=F, all.inside=F) {
  nx <- length(x)
  index <- integer(nx)
  .C('find_interv_vec', xt=as.double(vec), n=length(vec),
    x=as.double(x), nx=nx, as.logical(rightmost.closed),
    as.logical(all.inside), index, DUP = FALSE, NAOK=T,
    PACKAGE='base')
  index
}

> system.time(for (i in 1:10000) findInterval(a - 0.5, sort(b)))
   user  system elapsed 
   1.22    0.00    1.22 
> system.time(for (i in 1:10000) sapply(a, function(x)sum(b<x)))
   user  system elapsed 
   0.79    0.00    0.78 
> system.time(for (i in 1:10000) rowSums(outer(a, b, ">")))
   user  system elapsed 
   0.72    0.00    0.72 
> system.time(for (i in 1:10000) findInterval(a - 0.5, b[order(b)]))
   user  system elapsed 
   0.42    0.00    0.42 
> system.time(for (i in 1:10000) findInterval2(a - 0.5, b[order2(b)]))
   user  system elapsed 
   0.16    0.00    0.15 

定义findInterval2and的复杂性order2可能只有在您有大量 N 相当小的迭代时才有保证。

还有较大 N 的时间：

> a = rep(a, 100)
> b = rep(b, 100)
> system.time(for (i in 1:100) findInterval(a - 0.5, sort(b)))
   user  system elapsed 
   0.01    0.00    0.02 
> system.time(for (i in 1:100) sapply(a, function(x)sum(b<x)))
   user  system elapsed 
   0.67    0.00    0.68 
> system.time(for (i in 1:100) rowSums(outer(a, b, ">")))
   user  system elapsed 
   3.67    0.26    3.94 
> system.time(for (i in 1:100) findInterval(a - 0.5, b[order(b)]))
   user  system elapsed 
      0       0       0 
> system.time(for (i in 1:100) findInterval2(a - 0.5, b[order2(b)]))
   user  system elapsed 
      0       0       0 

Answer 2

一种选择是将outer()二元运算符函数>应用于aand b：

> outer(a, b, ">")
      [,1]  [,2]  [,3]  [,4]  [,5]
[1,] FALSE FALSE FALSE FALSE FALSE
[2,] FALSE FALSE  TRUE FALSE  TRUE
[3,]  TRUE FALSE  TRUE  TRUE  TRUE
[4,]  TRUE FALSE  TRUE  TRUE  TRUE
[5,]  TRUE  TRUE  TRUE  TRUE  TRUE
[6,]  TRUE  TRUE  TRUE  TRUE  TRUE
[7,] FALSE FALSE  TRUE FALSE  TRUE

Q 的答案由上面结果的行总和给出：

> rowSums(outer(a, b, ">"))
[1] 0 2 4 4 5 5 2

对于此示例数据集，此解决方案稍快一些findIntervals()，但也快不了多少：

> system.time(replicate(1000, findInterval(a - 0.5, sort(b))))
   user  system elapsed 
  0.131   0.000   0.132 
> system.time(replicate(1000, rowSums(outer(a, b, ">"))))
   user  system elapsed 
  0.078   0.000   0.079

它也比该sapply()版本稍快一些，但速度稍快一些：

> system.time(replicate(1000, sapply(a, function(x)sum(b<x))))
   user  system elapsed 
  0.082   0.000   0.082

@Charles 指出，示例中的大部分时间findInterval()都是由 所使用sort()，可以通过 来规避order()。完成此操作后，findInterval()解决方案比解决方案更快outer()：

> system.time(replicate(1000, findInterval(a - 0.5, b[order(b)])))
   user  system elapsed 
  0.049   0.000   0.049