Лучший способ выбрать случайные элементы из массива с минимальной разницей в R

Question

Я бы хотел случайным образом выбрать из массива определенное количество элементов таким образом, чтобы они всегда соблюдали ограничение в их взаимном расстоянии. Например, имея вектор a <- seq(1,1000), как я могу выбрать 20 элементов с минимальным расстоянием 15 между собой?

На данный момент я использую простую итерацию, для которой я отклоняю выбор всякий раз, когда он находится слишком близко к какому-либо элементу, но он громоздок и имеет тенденцию быть длинным, если число элементов для выбора велико. Есть ли передовая практика / функция для этого?

РЕДАКТИРОВАТЬ - Сводка ответов и анализа

До сих пор у меня было два рабочих ответа, которые я обернул в две конкретные функции.

# dash2 approach
# ---------------
rand_pick_min <- function(ar, min.dist, n.picks){
  stopifnot(is.numeric(min.dist), 
            is.numeric(n.picks), n.picks%%1 == 0)
  if(length(ar)/n.picks < min.dist) 
    stop('The number of picks exceeds the maximum number of divisions that the array allows which is: ', 
         floor(length(ar)/min.dist))
  picked <- array(NA, n.picks)
  copy <- ar
  for (i in 1:n.picks) {
    stopifnot(length(copy) > 0)  
    picked[i] <- sample(copy, 1)
    copy <- copy[ abs(copy - picked[i]) >= min.dist ]
  }
  return(picked)
}

# denis approach
# ---------------
rand_pick_min2 <- function(ar, min.dist, n.picks){
  require(Surrogate)
  stopifnot(is.numeric(min.dist), 
            is.numeric(n.picks), n.picks%%1 == 0)
  if(length(ar)/n.picks < min.dist) 
    stop('The number of picks exceeds the maximum number of divisions that the array allows which is: ', 
         floor(length(ar)/min.dist))
  lar <- length(ar)
  dist <- Surrogate::RandVec(a=min.dist, b=(lar-(n.picks)*min.dist), 
                             s=lar, n=(n.picks+1), m=1, Seed=sample(1:lar, size = 1))$RandVecOutput
  return(cumsum(round(dist))[1:n.picks])
}

Используя тот же пример, я запустил 3 теста. Во-первых, эффективная действительность минимального лимита

# Libs
require(ggplot2)
require(microbenchmark)

# Inputs
a <- seq(1, 1000)            # test vector
md <- 15                     # min distance
np <- 20                     # number of picks

# Run
dist_vec <- c(sapply(1:500, function(x) c(dist(rand_pick_min(a, md, np)))))   # sol 1
dist_vec2 <- c(sapply(1:500, function(x) c(dist(rand_pick_min2(a, md, np))))) # sol 2

# Tests - break the min
cat('Any distance breaking the min in sol 1?', any(dist_vec < md), '\n')  # FALSE
cat('Any distance breaking the min in sol 2?', any(dist_vec2 < md), '\n') # FALSE

Во-вторых, я проверил распределение полученных расстояний, получив первые два графика в порядке решения (sol1 [A] - это sol для dash2, а sol2 [B] - для denis ').

pa <- ggplot() + theme_classic() +
  geom_density(aes_string(x = dist_vec), fill = 'lightgreen') +
  geom_vline(aes_string(xintercept = mean(dist_vec)), col = 'darkred') + xlab('Distances')
pb <- ggplot() + theme_classic() +
  geom_density(aes_string(x = dist_vec2), fill = 'lightgreen') +
  geom_vline(aes_string(xintercept = mean(dist_vec)), col = 'darkred') + xlab('Distances')
print(pa)
print(pb)

Наконец, я вычислил время вычислений, необходимое для двух подходов, следующим образом и получив последнюю цифру.

comp_times <- microbenchmark::microbenchmark(
  'solution_1' = rand_pick_min(a, md, np),
  'solution_2' = rand_pick_min2(a, md, np),
  times = 500
)
ggplot2::autoplot(comp_times); ggsave('stckoverflow2.png')

Ободрившись результатами, я спрашиваю себя, должно ли быть ожидаемое распределение расстояний или это отклонение из-за примененных методов.

EDIT2 - Ответ на последний вопрос после комментария, сделанного Денисом

Используя намного больше процедур отбора проб (5000), я подготовил PDF-файл с результирующими позициями, и ваш подход действительно содержит некоторый артефакт, который заставляет ваше решение (B) отклоняться от того, которое мне было нужно. Тем не менее, было бы интересно иметь возможность обеспечить конкретное окончательное распределение позиций.

Answer 1

Если вы хотите избежать методов попадания и промаха, вам придется перевести вашу проблему в выборку расстояний с ограничениями на сумму ваших расстояний.

В основном, как я перевожу то, что вы хотите: ваши N выборочных позиций эквивалентны N + 1 расстоянию, в пределах от минимального расстояния до размера вашего вектора - N * mindist (случай, когда все ваши образцы упакованы вместе). Затем вам нужно ограничить сумму расстояний равной 1000 (размер вашего вектора).

В этом случае решение будет использовать пакет Surrogate :: RandVec из пакета Surrogate (см. Случайная выборка для получения точной суммы ), что позволяет производить выборку с фиксированной суммой.

library(Surrogate)
a <- seq(1,1000)
mind <- 15
N <- 20
dist <- Surrogate::RandVec(a=mind, b=(1000-(N)*mind), s=1000, n=(N+1), m=1, Seed=sample(1:1000, size = 1))$RandVecOutput
pos <- cumsum(round(dist))[1:20]
pos

> pos
 [1]  22  59  76 128 204 239 289 340 389 440 489 546 567 607 724 773 808 843 883 927

dist - выборка на расстоянии. Вы восстанавливаете свою позицию, делая сумму расстояний. Это дает вам pos, вектор ваших позиций индекса.

Преимущество состоит в том, что вы можете получить любое значение, и что ваша выборка должна быть случайной. Для части скорости, которую я не знаю, вам нужно сравнить с вашим методом для случая с большими данными.

Вот гистограмма 1000 попыток:

Answer 2

Я думаю, что лучшее решение, которое в некотором смысле гарантирует случайность (я не совсем уверен, какой смысл!) Может быть:

1. Выберите случайный элемент
2. Удалить все элементы, которые находятся слишком близко к этому элементу
3. Выберите другой элемент
4. Возврат к 2.

Итак:

min_dist <- 15
a <- seq(1, 1000)
picked <- integer(20)
copy <- a
for (i in 1:20) {
  stopifnot(length(copy) > 0)
  picked[i] <- sample(copy, 1)
  copy <- copy[ abs(copy - picked[i]) >= min_dist ]
}

Будет ли это быстрее, чем выборка и отбраковка, может зависеть от характеристик исходного вектора. Кроме того, как вы можете видеть, вы не гарантированно сможете получить все нужные вам элементы, хотя в вашем конкретном случае проблем не будет, потому что 19 интервалов шириной 30 никогда не смогут охватить весь seq(1, 1000).