2008-11-22 8 views
15

Próbuję sprawdzić prawdopodobieństwo, że dane klastrowanie danych wystąpiły przypadkowo. Mocnym sposobem na to jest symulacja Monte Carlo, w której powiązania między danymi i grupami są losowo przydzielane dużą liczbę razy (np. 10 000), a metryka grupowania jest używana do porównania rzeczywistych danych z symulacjami w celu określenia wartość.Algorytm do pobierania próbek bez wymiany?

Mam większość tego działającego, ze wskaźnikami mapującymi zgrupowanie na elementy danych, więc planuję losowe przypisywanie wskaźników do danych. PYTANIE: jaki jest szybki sposób próbkowania bez zastąpienia, tak aby każdy wskaźnik był losowo ponownie przypisywany w zestawach danych replik?

na przykład (dane te są tylko uproszczoną przykład):

danych (n = 12 wartości) - Grupa A: 0,1, 0,2, 0,4/Grupa B: 0,5, 0,6, 0,8/Grupa C : 0,4, 0,5/Grupa D: 0,2, 0,2, 0,3, 0,5

Dla każdego replikacji zestaw danych, ja te same rozmiary klastra (A = 3, B = 3, c = 2, d = 4) i wartości danych, ale ponownie przypisze wartości do klastrów.

Aby to zrobić, mogę wygenerować liczby losowe z zakresu 1-12, przypisać pierwszy element grupy A, a następnie wygenerować liczby losowe z zakresu 1-11 i przypisać drugi element w grupie A, i tak dalej . Zmiana przypisania wskaźnika jest szybka i wstępnie przydzielę wszystkie struktury danych, ale próbkowanie bez zastąpienia wydaje się być problemem, który mógł zostać rozwiązany wiele razy wcześniej.

Preferowana logika lub pseudokod.

Odpowiedz

5

Zobacz moją odpowiedź na to pytanie Unique (non-repeating) random numbers in O(1)?. Ta sama logika powinna osiągnąć to, czego szukasz.

+0

Świetnie! Niestety nie widziałem tej odpowiedzi, gdy szukałem SO (do pobierania próbek bez wymiany, statystyk, algorytmów itp.). Może to będzie służyć jako meta-pytanie, by poprowadzić ludzi takich jak ja do twojej oryginalnej odpowiedzi. Twoje zdrowie! – Argalatyr

35

Oto kod do próbkowania bez zastąpienia na podstawie Algorithm 3.4.2S z książki Knuth'a Seminumeric Algorithms.

void SampleWithoutReplacement 
(
    int populationSize, // size of set sampling from 
    int sampleSize,  // size of each sample 
    vector<int> & samples // output, zero-offset indicies to selected items 
) 
{ 
    // Use Knuth's variable names 
    int& n = sampleSize; 
    int& N = populationSize; 

    int t = 0; // total input records dealt with 
    int m = 0; // number of items selected so far 
    double u; 

    while (m < n) 
    { 
     u = GetUniform(); // call a uniform(0,1) random number generator 

     if ((N - t)*u >= n - m) 
     { 
      t++; 
     } 
     else 
     { 
      samples[m] = t; 
      t++; m++; 
     } 
    } 
} 

jest bardziej skuteczne, ale bardziej złożone metody Jeffrey Scott Vitter w "An Efficient algorytmu sekwencyjnego losowego wyboru, ACM Transactions on" Mathematical Software, 13 (1), marzec 1987, 58-67.

+0

Nie mam tej książki (jeszcze) i miałem problemy z udowodnieniem poprawności algorytmu dla siebie. Zaimplementowałem go w Javie i sprawdziłem, czy elementy populacji są próbkowane z jednolitym prawdopodobieństwem. Wyniki są przekonujące. Zobacz to [gist] (https://gist.github.com/10864989) – Alban

+0

Bezkrytyczna implementacja metody D Vittera w Mathematica jest o kilka rzędów wielkości szybsza niż wbudowany algorytm. Opiszę to tutaj: http://tinyurl.com/lbldlpq –

+4

@Alban - Możemy zobaczyć problem pobierania próbek n elementów z populacji N, biorąc pod uwagę pierwszy element.Istnieje prawdopodobieństwo (n/N), że ten element jest zawarty: jeśli tak, problem ogranicza się do elementów próbkowania (n-1) z pozostałych (N-1); jeśli nie, problem ogranicza się do próbkowania (n) elementów z pozostałych (N-1). Niektóre transformacje zmienne pokażą, że jest to istotą algorytmu Knutha (poprzez inkrementację t). –

0

Inny algorytm pobierania próbek bez zastąpienia jest opisany here.

Jest podobny do opisanego przez Johna D. Cooka w jego odpowiedzi, a także od Knutha, ale ma inną hipotezę: Rozmiar populacji jest nieznany, ale próbka może zmieścić się w pamięci. Ten nazywa się "algorytmem Knutha S".

Cytując artykuł rosettacode:

  1. wybierz pierwszą pozycje n jako próbki, jak stają się one dostępne;
  2. Dla i-tego elementu, w którym i> n, istnieje losowa szansa na utrzymanie go w stanie nienaruszonym. Jeśli nie uda się tej szansy, próbka pozostanie taka sama. Jeśli nie, to niech losowo (1/n) zastąpi jedną z wcześniej wybranych pozycji próbki.
  3. Powtórz 2 dla kolejnych pozycji.
+1

Rosettacode ma złą nazwę dla algorytmu: powinien to być "Algorytm R" lub "Pobieranie próbek z rezerwuaru". "Algorytm S" (inaczej "Wybór metody próbkowania") wymaga wcześniejszego określenia wielkości populacji. Oba algorytmy są opisane w TAOCP - Tom 2 - § 4.3.2 – manlio

7

A C++ kod działa w oparciu o answer by John D. Cook.

#include <random> 
#include <vector> 

double GetUniform() 
{ 
    static std::default_random_engine re; 
    static std::uniform_real_distribution<double> Dist(0,1); 
    return Dist(re); 
} 

// John D. Cook, https://stackoverflow.com/a/311716/15485 
void SampleWithoutReplacement 
(
    int populationSize, // size of set sampling from 
    int sampleSize,  // size of each sample 
    std::vector<int> & samples // output, zero-offset indicies to selected items 
) 
{ 
    // Use Knuth's variable names 
    int& n = sampleSize; 
    int& N = populationSize; 

    int t = 0; // total input records dealt with 
    int m = 0; // number of items selected so far 
    double u; 

    while (m < n) 
    { 
     u = GetUniform(); // call a uniform(0,1) random number generator 

     if ((N - t)*u >= n - m) 
     { 
      t++; 
     } 
     else 
     { 
      samples[m] = t; 
      t++; m++; 
     } 
    } 
} 

#include <iostream> 
int main(int,char**) 
{ 
    const size_t sz = 10; 
    std::vector<int> samples(sz); 
    SampleWithoutReplacement(10*sz,sz,samples); 
    for (size_t i = 0; i < sz; i++) { 
    std::cout << samples[i] << "\t"; 
    } 

    return 0; 
} 
2

Zainspirowany przez @John D. Cook's answer, napisałem implementację w Nim. Na początku miałem trudności ze zrozumieniem, jak to działa, więc obszernie skomentowałem również przykład. Może to pomaga zrozumieć ten pomysł. Zmieniłem też nieznacznie nazwy zmiennych.

iterator uniqueRandomValuesBelow*(N, M: int) = 
    ## Returns a total of M unique random values i with 0 <= i < N 
    ## These indices can be used to construct e.g. a random sample without replacement 
    assert(M <= N) 

    var t = 0 # total input records dealt with 
    var m = 0 # number of items selected so far 

    while (m < M): 
    let u = random(1.0) # call a uniform(0,1) random number generator 

    # meaning of the following terms: 
    # (N - t) is the total number of remaining draws left (initially just N) 
    # (M - m) is the number how many of these remaining draw must be positive (initially just M) 
    # => Probability for next draw = (M-m)/(N-t) 
    # i.e.: (required positive draws left)/(total draw left) 
    # 
    # This is implemented by the inequality expression below: 
    # - the larger (M-m), the larger the probability of a positive draw 
    # - for (N-t) == (M-m), the term on the left is always smaller => we will draw 100% 
    # - for (N-t) >> (M-m), we must get a very small u 
    # 
    # example: (N-t) = 7, (M-m) = 5 
    # => we draw the next with prob 5/7 
    # lets assume the draw fails 
    # => t += 1 => (N-t) = 6 
    # => we draw the next with prob 5/6 
    # lets assume the draw succeeds 
    # => t += 1, m += 1 => (N-t) = 5, (M-m) = 4 
    # => we draw the next with prob 4/5 
    # lets assume the draw fails 
    # => t += 1 => (N-t) = 4 
    # => we draw the next with prob 4/4, i.e., 
    # we will draw with certainty from now on 
    # (in the next steps we get prob 3/3, 2/2, ...) 
    if (N - t)*u >= (M - m).toFloat: # this is essentially a draw with P = (M-m)/(N-t) 
     # no draw -- happens mainly for (N-t) >> (M-m) and/or high u 
     t += 1 
    else: 
     # draw t -- happens when (M-m) gets large and/or low u 
     yield t # this is where we output an index, can be used to sample 
     t += 1 
     m += 1 

# example use 
for i in uniqueRandomValuesBelow(100, 5): 
    echo i 
1

Gdy rozmiar populacji jest znacznie większy niż rozmiar próbki wyżej algorytmy się niewydajne, ponieważ nie mają złożoność O (n) n jest wielkość populacji.

Kiedy byłem studentem pisałem kilka algorytmów jednolitego próbkowania bez wymiany, które mają średnią złożoność O (s dziennika s), gdzie s jest wielkość próbki. Oto kod algorytmu drzewo binarne, przy średniej złożoności O (s log s) w R:

# The Tree growing algorithm for uniform sampling without replacement 
# by Pavel Ruzankin 
quicksample = function (n,size) 
# n - the number of items to choose from 
# size - the sample size 
{ 
    s=as.integer(size) 
    if (s>n) { 
    stop("Sample size is greater than the number of items to choose from") 
    } 
    # upv=integer(s) #level up edge is pointing to 
    leftv=integer(s) #left edge is poiting to; must be filled with zeros 
    rightv=integer(s) #right edge is pointig to; must be filled with zeros 
    samp=integer(s) #the sample 
    ordn=integer(s) #relative ordinal number 

    ordn[1L]=1L #initial value for the root vertex 
    samp[1L]=sample(n,1L) 
    if (s > 1L) for (j in 2L:s) { 
    curn=sample(n-j+1L,1L) #current number sampled 
    curordn=0L #currend ordinal number 
    v=1L #current vertice 
    from=1L #how have come here: 0 - by left edge, 1 - by right edge 
    repeat { 
     curordn=curordn+ordn[v] 
     if (curn+curordn>samp[v]) { #going down by the right edge 
     if (from == 0L) { 
      ordn[v]=ordn[v]-1L 
     } 
     if (rightv[v]!=0L) { 
      v=rightv[v] 
      from=1L 
     } else { #creating a new vertex 
      samp[j]=curn+curordn 
      ordn[j]=1L 
      # upv[j]=v 
      rightv[v]=j 
      break 
     } 
     } else { #going down by the left edge 
     if (from==1L) { 
      ordn[v]=ordn[v]+1L 
     } 
     if (leftv[v]!=0L) { 
      v=leftv[v] 
      from=0L 
     } else { #creating a new vertex 
      samp[j]=curn+curordn-1L 
      ordn[j]=-1L 
      # upv[j]=v 
      leftv[v]=j 
      break 
     } 
     } 
    } 
    } 
    return(samp) 
} 

złożoności tego algorytmu opisane w: Rouzankin P. S .; Voytishek, A. V. O kosztach algorytmów doboru losowego. Monte Carlo Methods Appl. 5 (1999), nie. 1, 39-54. http://dx.doi.org/10.1515/mcma.1999.5.1.39

Jeśli uznasz, że algorytm jest przydatny, prosimy o odniesienie.

Zobacz także: P. Gupta, G. P. Bhattacharjee. (1984) Skuteczny algorytm losowego pobierania próbek bez wymiany. International Journal of Computer Mathematics 16: 4, strony 201-209. DOI: 10.1080/00207168408803438

Teuhola, J. and Nevalainen, O. 1982. Dwa wydajne algorytmy losowego pobierania próbek bez wymiany./IJCM /, 11 (2): 127-140. DOI: 10,1080/00207168208803304

W ostatniej pracy autorzy wykorzystać tabele mieszania i twierdzenia, że ​​algorytmy O (s) złożoności. Istnieje jeszcze jeden algorytm szybkiego mieszania tablic, który wkrótce zostanie zaimplementowany w pqR (dość szybki R): https://stat.ethz.ch/pipermail/r-devel/2017-October/075012.html