Como reamostrar em R sem repetir permutações?

12

Em R, se eu definir.seed () e usar a função de amostra para randomizar uma lista, posso garantir que não gerarei a mesma permutação?

ie ...

set.seed(25)
limit <- 3
myindex <- seq(0,limit)
for (x in seq(1,factorial(limit))) {
    permutations <- sample(myindex)
    print(permutations)
}

Isso produz

[1] 1 2 0 3
[1] 0 2 1 3
[1] 0 3 2 1
[1] 3 1 2 0
[1] 2 3 0 1
[1] 0 1 3 2

todas as permutações impressas serão permutações únicas? Ou existe alguma chance, com base na maneira como isso é implementado, de obter algumas repetições?

Quero poder fazer isso sem repetições, garantido. Como eu faria isso?

(Também quero evitar ter que usar uma função como permn (), que tem um método muito mecanicista para gerar todas as permutações - não parece aleatório.)

Além disso, nota de rodapé --- parece que esse problema é O ((n!)!), Se não me engano.

Mittenchops
fonte
Por padrão, o argumento 'substituir' de 'amostra' é definido como FALSE.
Ocram
Obrigado ocram, mas está funcionando dentro de uma amostra específica. Portanto, isso garante que 0,1,2 e 3 não se repitam em um empate (portanto, não posso desenhar 0,1,2,2), mas não sei se isso garante que a segunda amostra, Não consigo desenhar a mesma sequência de 0123 novamente. É isso que estou pensando em termos de implementação, se definir a semente tem algum efeito nessa repetição.
Mittenchops 8/03/12
Sim, é isso que eu finalmente entendi lendo as respostas ;-)
Ocram
1
Se limitexceder 12, você provavelmente ficará sem memória RAM quando R tentar alocar espaço para seq(1,factorial(limit)). (! 12 requer cerca de 2 GB, de modo 13 terá cerca de 25 GB, 14 cerca de 350 GB, etc.!)
whuber
2
Existe uma solução rápida, compacta e elegante para gerar seqüências aleatórias de todas as permutações de 1: n, desde que você possa armazenar confortavelmente n! números inteiros no intervalo 0: (n!). Combina a representação da tabela de inversão de uma permutação com a representação da base fatorial de números.
whuber

Respostas:

9

A questão tem muitas interpretações válidas. Os comentários - especialmente aquele que indica permutações de 15 ou mais elementos são necessários (15! = 1307674368000 está ficando grande) - sugerem que o que se quer é uma amostra aleatória relativamente pequena , sem substituição, de todos os n! = n * (n-1) (n-2) ... * 2 * 1 permutações de 1: n. Se isso for verdade, existem soluções (um tanto) eficientes.

A função a seguir rperm, aceita dois argumentos n(o tamanho das permutações para amostrar) e m(o número de permutações de tamanho n para desenhar). Se m se aproximar ou exceder n !, a função levará muito tempo e retornará muitos valores de NA: ela deve ser usada quando n for relativamente grande (digamos, 8 ou mais) e m for muito menor que n !. Ele funciona armazenando em cache uma representação em cadeia das permutações encontradas até o momento e gerando novas permutações (aleatoriamente) até que uma nova seja encontrada. Ele explora a capacidade de indexação de lista associativa de R para pesquisar rapidamente a lista de permutações encontradas anteriormente.

rperm <- function(m, size=2) { # Obtain m unique permutations of 1:size

    # Function to obtain a new permutation.
    newperm <- function() {
        count <- 0                # Protects against infinite loops
        repeat {
            # Generate a permutation and check against previous ones.
            p <- sample(1:size)
            hash.p <- paste(p, collapse="")
            if (is.null(cache[[hash.p]])) break

            # Prepare to try again.
            count <- count+1
            if (count > 1000) {   # 1000 is arbitrary; adjust to taste
                p <- NA           # NA indicates a new permutation wasn't found
                hash.p <- ""
                break
            }
        }
        cache[[hash.p]] <<- TRUE  # Update the list of permutations found
        p                         # Return this (new) permutation
    }

    # Obtain m unique permutations.
    cache <- list()
    replicate(m, newperm())  
} # Returns a `size` by `m` matrix; each column is a permutation of 1:size.

A natureza de replicateé retornar as permutações como vetores de coluna ; por exemplo , o seguinte reproduz um exemplo na pergunta original, transposta :

> set.seed(17)
> rperm(6, size=4)
     [,1] [,2] [,3] [,4] [,5] [,6]
[1,]    1    2    4    4    3    4
[2,]    3    4    1    3    1    2
[3,]    4    1    3    2    2    3
[4,]    2    3    2    1    4    1

Os tempos são excelentes para valores pequenos a moderados de m, até cerca de 10.000, mas degradam-se para problemas maiores. Por exemplo, uma amostra de m = 10.000 permutações de n = 1000 elementos (uma matriz de 10 milhões de valores) foi obtida em 10 segundos; uma amostra de m = 20.000 permutações de n = 20 elementos exigiu 11 segundos, mesmo que a saída (uma matriz de 400.000 entradas) fosse muito menor; e a amostra de computação de m = 100.000 permutações de n = 20 elementos foi abortada após 260 segundos (não tive paciência para aguardar a conclusão). Esse problema de escala parece estar relacionado a ineficiências de escala no endereçamento associativo de R. Pode-se contornar isso gerando amostras em grupos de, digamos, mais ou menos 1000, combinando essas amostras em uma amostra grande e removendo duplicatas.

Editar

Podemos obter um desempenho assintótico quase linear dividindo o cache em uma hierarquia de dois caches, para que R nunca precise pesquisar em uma lista grande. Conceitualmente (embora não seja implementado), crie uma matriz indexada pelos primeiros elementos de uma permutação. As entradas nessa matriz são listas de todas as permutações que compartilham esses primeiros elementos. Para verificar se uma permutação foi vista, use seus primeiros elementos para encontrar sua entrada no cache e, em seguida, procure essa permutação nessa entrada. Podemos escolher para equilibrar os tamanhos esperados de todas as listas. A implementação real não usa umk k k kkkkkkarray -fold, que seria difícil de programar em generalidade suficiente, mas usa outra lista.

Aqui estão alguns tempos decorridos em segundos para uma variedade de tamanhos de permutação e número de permutações distintas solicitadas:

 Number Size=10 Size=15 Size=1000 size=10000 size=100000
     10    0.00    0.00      0.02       0.08        1.03
    100    0.01    0.01      0.07       0.64        8.36
   1000    0.08    0.09      0.68       6.38
  10000    0.83    0.87      7.04      65.74
 100000   11.77   10.51     69.33
1000000  195.5   125.5

(A aceleração aparentemente anômala de size = 10 para size = 15 é porque o primeiro nível do cache é maior para size = 15, reduzindo o número médio de entradas nas listas de segundo nível, acelerando a pesquisa associativa de R. custo em RAM, a execução pode ser feita mais rapidamente aumentando o tamanho do cache de nível superior, aumentando apenas k.head1 (o que multiplica o tamanho de nível superior por 10) acelerando rperm(100000, size=10)de 11,77 segundos para 8,72 segundos, por exemplo. cache 10 vezes maior, mas não obteve ganho apreciável, com clock de 8,51 segundos.)

Exceto no caso de 1.000.000 de permutações únicas de 10 elementos (uma parcela substancial de todos os 10! = Cerca de 3,63 milhões de permutações), praticamente nenhuma colisão foi detectada. Nesse caso excepcional, houve 169.301 colisões, mas nenhuma falha completa (um milhão de permutações únicas foram de fato obtidas).

Observe que, com tamanhos de permutação grandes (maiores que 20 ou mais), a chance de obter duas permutações idênticas, mesmo em uma amostra de até 1.000.000.000, é muito pequena. Assim, esta solução é aplicável principalmente em situações em que (a) um grande número de permutações únicas de (b) entre e ou então elementos estão a ser gerada, mas mesmo assim, (c), substancialmente menos do que todos ospermutações são necessárias.n = 15 n !n=5n=15n!

Código de trabalho a seguir.

rperm <- function(m, size=2) { # Obtain m unique permutations of 1:size
    max.failures <- 10

    # Function to index into the upper-level cache.
    prefix <- function(p, k) {    # p is a permutation, k is the prefix size
        sum((p[1:k] - 1) * (size ^ ((1:k)-1))) + 1
    } # Returns a value from 1 through size^k

    # Function to obtain a new permutation.
    newperm <- function() {
        # References cache, k.head, and failures in parent context.
        # Modifies cache and failures.        

        count <- 0                # Protects against infinite loops
        repeat {
            # Generate a permutation and check against previous ones.
            p <- sample(1:size)
            k <- prefix(p, k.head)
            ip <- cache[[k]]
            hash.p <- paste(tail(p,-k.head), collapse="")
            if (is.null(ip[[hash.p]])) break

            # Prepare to try again.
            n.failures <<- n.failures + 1
            count <- count+1
            if (count > max.failures) {  
                p <- NA           # NA indicates a new permutation wasn't found
                hash.p <- ""
                break
            }
        }
        if (count <= max.failures) {
            ip[[hash.p]] <- TRUE      # Update the list of permutations found
            cache[[k]] <<- ip
        }
        p                         # Return this (new) permutation
    }

    # Initialize the cache.
    k.head <- min(size-1, max(1, floor(log(m / log(m)) / log(size))))
    cache <- as.list(1:(size^k.head))
    for (i in 1:(size^k.head)) cache[[i]] <- list()

    # Count failures (for benchmarking and error checking).
    n.failures <- 0

    # Obtain (up to) m unique permutations.
    s <- replicate(m, newperm())
    s[is.na(s)] <- NULL
    list(failures=n.failures, sample=matrix(unlist(s), ncol=size))
} # Returns an m by size matrix; each row is a permutation of 1:size.
whuber
fonte
Isso está próximo, mas noto que recebo alguns erros, como 1, 2 e 4, mas acho que entendi o que você quer dizer e deve poder trabalhar com ele. Obrigado! > rperm(6,3) $failures [1] 9 $sample [,1] [,2] [,3] [1,] 3 1 3 [2,] 2 2 1 [3,] 1 3 2 [4,] 1 2 2 [5,] 3 3 1 [6,] 2 1 3
Mittenchops
3

Usar uniqueda maneira correta deve fazer o truque:

set.seed(2)
limit <- 3
myindex <- seq(0,limit)

endDim<-factorial(limit)
permutations<-sample(myindex)

while(is.null(dim(unique(permutations))) || dim(unique(permutations))[1]!=endDim) {
    permutations <- rbind(permutations,sample(myindex))
}
# Resulting permutations:
unique(permutations)

# Compare to
set.seed(2)
permutations<-sample(myindex)
for(i in 1:endDim)
{
permutations<-rbind(permutations,sample(myindex))
}
permutations
# which contains the same permutation twice
MånsT
fonte
Desculpe por não explicar o código corretamente. Estou com pressa agora, mas fico feliz em responder a quaisquer perguntas que você tiver mais tarde. Além disso, eu não tenho nenhuma idéia sobre a velocidade do código acima ...
MånsT
1
Funcionalizei o que você me deu dessa maneira: `myperm <- function (limit) {myindex <- seq (0, limit) endDim <- fatorial (limit) permutações <- sample (myindex) while (is.null (dim (unique (permutações))) || dim (único (permutações)) [1]! = endDim) {permutações <- rbind (permutações, amostra (myindex))} retorno (único (permutações))} 'Funciona, mas enquanto eu pode fazer limit = 6, limit = 7 faz com que meu computador superaqueça. = PI acho que ainda deve haver uma maneira de subamostra isso ...
Mittenchops
@ Mittenchops, por que você diz que precisamos usar o exclusivo para reamostrar em R sem repetir permutações? Obrigado.
Frank
2

Vou desviar um pouco a sua primeira pergunta e sugerir que, se você está lidando com vetores relativamente curtos, você pode simplesmente gerar todas as permutações usando permne ordená-las aleatoriamente as que usam sample:

x <- combinat:::permn(1:3)
> x[sample(factorial(3),factorial(3),replace = FALSE)]
[[1]]
[1] 1 2 3

[[2]]
[1] 3 2 1

[[3]]
[1] 3 1 2

[[4]]
[1] 2 1 3

[[5]]
[1] 2 3 1

[[6]]
[1] 1 3 2
joran
fonte
Eu gosto muito disso, e tenho certeza de que é o pensamento certo. Mas meu problema me faz usar uma sequência que vai até 10. Permn () foi significativamente mais lento entre fatorial (7) e fatorial (8), então eu acho que 9 e 10 serão proibitivamente grandes.
Mittenchops
@Mittenchops True, mas ainda é possível que você realmente só precise calculá-las uma vez, certo? Salve-os em um arquivo e, em seguida, carregue-os quando precisar deles e faça "amostra" da lista predefinida. Então você pode fazer o cálculo lento permn(10)ou o que quer que seja apenas uma vez.
joran
Certo, mas se estou armazenando todas as permutações em algum lugar, mesmo isso se decompõe em fatorial (15) --- simplesmente muito espaço para armazenar. É por isso que estou me perguntando se definir a semente me permitirá coletar permutações coletivamente - e, se não, se existe um algoritmo para isso.
Mittenchops
@Mittenchops Definir uma semente não influenciará o desempenho, apenas garante o mesmo início toda vez que você fizer uma chamada para o PRNG.
Roman Luštrik 15/03/12
1
@Mitten Consulte a ajuda para set.seed: descreve como salvar o estado do RNG e restaurá-lo mais tarde.
whuber