获取匹配条件的第一个值（循环太慢）

Question

我有很多与此类似的矩阵，但是有数千行：

r <- 10
c <- 2
set.seed(333)

m1 <- matrix(runif(r*c)+1, r, c)

> m1
          [,1]     [,2]
 [1,] 1.467001 1.393902
 [2,] 1.084598 1.474218
 [3,] 1.973485 1.891222
 [4,] 1.571306 1.665011
 [5,] 1.020119 1.736832
 [6,] 1.723557 1.911469
 [7,] 1.609394 1.637850
 [8,] 1.306719 1.864651
 [9,] 1.063510 1.287575
[10,] 1.305353 1.129959

我有一个循环，告诉我，对于第一列的每个值，第二列中第一个值的索引是多少，该索引高10％像这样：

result <- 1:nrow(m1)

for (i in 1:nrow(m1)){
    result[i] <- which(m1[,2]>(1.1*m1[,1][i]))[1]
}
> result
 [1]  3  1 NA  3  1  6  3  2  1  2

我的矩阵太多了，要花几个小时，在对我的代码进行性能分析后，迄今为止最大的耗时任务是此循环。根据您的说法，最快的方法是什么？

例如，r = 30000：

start_time <- Sys.time()

for (i in 1:nrow(m1)){
    result[i] <- which(m1[,2]>(1.1*m1[,1][i]))[1]
}

end_time <- Sys.time()
a <- end_time - start_time

> a
Time difference of 11.25815 secs

谢谢您的帮助！

Answer 1

您可以在此处使用一些快捷方式。您正在寻找第2列中的第一个值大于其他值的第一个值。这意味着，永远不要去看比我们先前在第2栏中看到的值低的值。

在具有10行的示例中，将如下所示：

> cummax(m1[, 2])
 [1] 1.393902 1.474218 1.891222 1.891222 1.891222 1.911469 1.911469 1.911469 1.911469 1.911469
> which(cummax(m1[, 2]) == m1[, 2])
[1] 1 2 3 6

如您所见，这些是结果向量中的唯一值。

可以进行的第二个优化是对第一列进行排序。如果您首先开始寻找最低值并逐步提高，则不必每次都浏览第二列。如果与左行不再匹配，则只需要转到下一行。

这确实承担了对矩阵进行排序的成本，但是之后可以使用一次遍历两列来找到结果。

dostuff <- function(m1){
  orderColumn1 <- order(m1[, 1])

  plus.10 <- m1[, 1] * 1.1

  results <- rep(NA, length(plus.10))

  IndexColumn1 <- 1
  IndexColumn2 <- 1
  row2CurrentMax <- 0
  while(IndexColumn2 <= nrow(m1)){
    row2Current <- m1[IndexColumn2, 2]
    if(row2Current > row2CurrentMax){
      row2CurrentMax <- row2Current
      while(TRUE){
        row1Current <- plus.10[orderColumn1[IndexColumn1]]
        if(row1Current <= row2CurrentMax){
          results[orderColumn1[IndexColumn1]] <- IndexColumn2
          IndexColumn1 <- IndexColumn1 + 1
        } else {
          break
        }
      }
    }
    IndexColumn2 <- IndexColumn2 + 1
  }
  results
}

30000行：

> result <- dostuff(m1)
> end_time <- Sys.time()
> a <- end_time - start_time
> a
Time difference of 0.0600059 secs

Answer 2

我不认为这是最快的方法，但是会比使用当前的for循环方法更快。

plus.10 <- m1[, 1] * 1.1
m2 <- m1[,2]
result <- sapply( plus.10, function(x) which.min(m2 < x))
result[plus.10 > max(m2) ] <- NA

result
[1]  3  1 NA  3  1  6  3  2  1  2

编辑：根据Ronak的要求，microbenchmark到目前为止在10000行上提出的解决方案的结果：

Unit: milliseconds
   expr        min        lq       mean      median          uq         max neval   cld
     h1 335.342689 337.35915 361.320461  341.804840  347.856556  516.230972    25  b   
 sindri 672.587291 688.78673 758.445467  713.240778  811.298608 1049.109844    25    d 
     op 865.567412 884.99514 993.066179 1006.694036 1026.434344 1424.755409    25     e
   loco 675.809092 682.98591 731.256313  693.672064  807.007358  821.893865    25    d 
 dmitry 420.869493 427.56492 454.439806  433.656519  438.367480  607.030825    25   c  
    jad   4.369628   4.41044   4.735393    4.503657    4.556527    7.488471    25 a