Самый быстрый способ построения последовательности `c(1:1, 1:2, …, 1: n)`

#r #performance

#r #Производительность

Вопрос:

Для данного положительного n целого числа я хочу знать самый быстрый базовый алгоритм R (not Rcpp ) для построения целочисленного вектора c(1:1, 1:2, ..., 1:n) , который имеет длину n*(n 1)/2 . Есть бонусные баллы за быстрые и эффективные с точки зрения памяти алгоритмы, поскольку я в конечном итоге хочу избежать выделения вектора длины n*n .

Я знаю как минимум о двух подходах:

  • unlist(lapply(seq_len(n), seq_len), FALSE, FALSE)
  • {J <- .row(c(n, n)); J[upper.tri(J, TRUE)]}

последнее особенно неэффективно, поскольку оно выделяет два целых вектора длины n*n .

Обратите внимание, что если мы присвоим значение .col(c(n, n)) J выше, то получим последовательность 1 2 2 3 3 3 4 4 4 4 ... . Эта последовательность может быть построена быстро и эффективно с {i <- seq_len(n); rep.int(i, i)} помощью .

Мне интересно, существует ли в этом случае такой же быстрый (или более быстрый) алгоритм .row(c(n, n)) , или если unlist lapply является оптимальным с точки зрения базового R.

FWIW, вот эталон из трех процедур, о которых я упоминал до сих пор:

 ## Seemingly optimal for 1 2 2 3 3 3 4 4 4 4 ...
f0 <- function(n) {i <- seq_len(n); rep.int(i, i)}
## Candidates for 1 1 2 1 2 3 1 2 3 4 ... (the sequence I actually want)
f1 <- function(n) unlist(lapply(seq_len(n), seq_len), FALSE, FALSE)
f2 <- function(n) {J <- .row(c(n, n)); J[upper.tri(J, TRUE)]}

n <- 1000L
microbenchmark::microbenchmark(f0(n), f1(n), f2(n), times = 10000L)
 
 Unit: milliseconds
  expr      min       lq     mean   median       uq      max neval
 f0(n) 1.711873 1.797891 2.112043 1.810273 1.836636 14.96644 10000
 f1(n) 1.986737 2.108630 2.472612 2.148195 2.214369 15.16282 10000
 f2(n) 3.785981 4.624821 5.551115 5.051405 5.861954 17.28740 10000
 

(Я знаю, что f1 это довольно близко к f0 этому, но есть ли что-то лучше, чем f1 ?)

Комментарии:

1. c(1, 1, 2, 1, 2, 3, 1, 2, 3, 4)

2. Я исправил заголовок…

3. Я думаю, вы могли бы попробовать Rcpp, если вам нужно быстрое и экономичное решение.

4. f0 не идентичен f1. f1 получает желаемый результат, но f0 уходит 1 2 2 3 3 3 ... . Тем не менее, я бы ожидал, что решение, подобное f0, будет самым быстрым, доступным в базе R, хотя бы потому, что оно переносит все циклы на C . «Медлительность», вероятно, возникает из-за принудительного выделения памяти при преобразовании последовательностей ALTREP в последовательности в памяти посредством конкатенации, в основном один из следующих форматов является компактным, а другой — нет : .Internal(inspect(1:10)); .Internal(inspect(f1(10))) . «Самым быстрым» было бы создать новый класс ALTREP на C / C для этой конкретной цели.

5. Я думаю sequences , что это самое сложное решение R уже, если вы не хотите использовать Rcpp

Ответ №1:

Я не уверен, о чем вы знаете, но если функция from base подходит, попробуйте sequence .

 f3 <- function(n) {sequence(1:n)}
 

Кажется, это почти в 2 ~ 3 раза быстрее, чем f0

Комментарии:

1. Спасибо — я не знал об этой функции. Похоже, это правильный базовый R-способ получить именно тот результат, который я хочу.


Ответ №2:

Я думаю sequence , это тот, который вам нужен (если вы не собираетесь использовать Rcpp для еще более быстрой версии)

 f1 <- function(n) unlist(lapply(seq_len(n), seq_len), FALSE, FALSE)
f2 <- function(n) {
  J <- .row(c(n, n))
  J[upper.tri(J, TRUE)]
}
f3 <- function(n) {
  v <- 1:n
  data.table::rowid(rep.int(v, v))
}
f4 <- function(n) sequence(1:n)

n <- 1000L
microbenchmark::microbenchmark(f1(n), f2(n), f3(n), f4(n), check = "identical")
 

Сравнительный анализ

 > microbenchmark::microbenchmark(f1(n), f2(n), f3(n), f4(n), check = "identical")
Unit: microseconds
  expr    min       lq      mean  median       uq     max neval
 f1(n) 3928.8  4144.50  5185.839  4227.5  4289.15 67457.1   100
 f2(n) 9490.3 10083.90 14415.777 12951.0 15080.50 78014.2   100
 f3(n) 8083.5  8572.10 12154.922  9063.0  9534.45 75408.7   100
 f4(n)  213.9   425.05   787.637   442.6   494.00  7844.4   100
 

Комментарии:

1. Я не думаю, что мой сеанс знает rowid .

2. @MikaelJagan Извините, что я забыл упомянуть, что это из data.table пакета

3. Если это так, удаленный ответ Парка был правильным, давайте проголосуем за восстановление.

4. @zx8754 да, это уже сделано

5. Я принял предыдущий ответ, который, как я полагаю, был невидим для меня, поскольку у меня нет привилегий. Спасибо за ваши комментарии / исправления.

Ответ №3:

Эти 2 также могут быть вариантами-

 n <- 5

unlist(purrr::map(seq(5), ~seq(.x)))
#>  [1] 1 1 2 1 2 3 1 2 3 4 1 2 3 4 5

unlist(mapply(FUN = function(.x) seq(.x), seq(n)))
#>  [1] 1 1 2 1 2 3 1 2 3 4 1 2 3 4 5
 

Создано 2021-12-10 пакетом reprex (версия v2.0.1)