Получить первое значение, которое соответствует условию (слишком медленный цикл) - программирование

Получить первое значение, которое соответствует условию (слишком медленный цикл)

У меня есть много матриц, похожих на это, но с тысячами строк:

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

Спасибо за помощь!

4b9b3361

ОТВЕТЫ

Ответ 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

Ответ 2

Я не думаю, что это самый быстрый способ, но он будет несколько быстрее, чем использование токового цикла.

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

Edit: В соответствии с просьбой Ронак, 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

Ответ 3

Вот попытка использования match() которая сокращает время по сравнению с примером r = 30000 в исходном сообщении примерно на 25%.

sapply(m1[, 1] * 1.1, function(x) match(TRUE, m1[, 2] > x))

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

Ответ 4

Лучший способ оптимизировать ваш код - использовать пакет data.table

Этот код дает вам> 2-кратное ускорение.

library(data.table);

setDTthreads(0);

r <- 30000;
c <- 2;
set.seed(333);

m1 <- matrix(runif(r*c)+1, r, c);
result1 <- rep(NA, nrow(m1));

start_time <- Sys.time();

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

#result1

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


start_time <- Sys.time()

tstDT <- data.table(m1);
#result2 <- tstDT[, sapply(V1, function(elem) { which(V2 > 1.1*elem)[1] })]
result2 <- tstDT[, sapply(V1, function(x) match(TRUE, V2 > 1.1*x) )]

#result2

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

Небольшой комментарий - я использую data.table, скомпилированный gcc с march = native и O3. Возможные O2 и march = core (как в стандартном пакете при установке) ускорения будут меньше, но...

Результат:

> library(data.table);
> 
> setDTthreads(0);
> 
> r <- 30000;
> c <- 2;
> set.seed(333);
> 
> m1 <- matrix(runif(r*c)+1, r, c);
> result1 <- rep(NA, nrow(m1));
> 
> start_time <- Sys.time();
> 
> for (i in 1:nrow(m1))
+ {
+     result1[i] <- which(m1[,2]>(1.1*m1[,1][i]))[1];
+ }
> 
> #result1
> 
> end_time <- Sys.time()
> a <- end_time - start_time
> a
Time difference of 8.738938 secs
> 
> 
> start_time <- Sys.time()
> 
> tstDT <- data.table(m1);
> #result2 <- tstDT[, sapply(V1, function(elem) { which(V2 > 1.1*elem)[1] })]
> result2 <- tstDT[, sapply(V1, function(x) match(TRUE, V2 > 1.1*x) )]
> 
> #result2
> 
> end_time <- Sys.time()
> a <- end_time - start_time
> a
Time difference of 3.582921 secs
> 
> 
> 
>

Ответ 5

Я предлагаю это:

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

m1 <- matrix(runif(r*c)+1, r, c)
x2 <-m1[, 2]



start_time <- Sys.time()

result <- lapply(m1[, 1], function(x) {
  min(which(m1[,2]>(1.1*x)))
})
end_time <- Sys.time()
a <- end_time - start_time
a


start_time <- Sys.time()

result <- lapply(m1[, 1], function(x) {
            min(which(x2>(1.1*x)))
})
end_time <- Sys.time()
a <- end_time - start_time
a

Первый: 8,6 с Второй: 6,4 с