Получить первое значение, которое соответствует условию (слишком медленный цикл)
У меня есть много матриц, похожих на это, но с тысячами строк:
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
Спасибо за помощь!
Ответы
Ответ 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 с