Generuj losowe liczby po rozkładzie w przedziale

17

Muszę wygenerować liczby losowe po rozkładzie normalnym w przedziale $(a,b)$ . (Pracuję w R.)

Wiem, że funkcja rnorm(n,mean,sd)wygeneruje losowe liczby po rozkładzie normalnym, ale jak ustawić limity interwałów w tym zakresie? Czy są do tego dostępne jakieś konkretne funkcje R?

— dvs
źródło

Dlaczego chcesz to zrobić? Jeśli jest ograniczony, to naprawdę nie może być normalny. Co próbujesz osiągnąć?

— gung - Przywróć Monikę

x <- rnorm(n, mean, sd); x <- x[x > lower.limit & x < upper.limit]

— Hugh

3

@Hugh, to świetnie ... o ile nie dbasz o liczbę losowych wartości.

— Glen_b

31

To brzmi jak chcesz symulować z obciętego rozkładu , aw twoim konkretnym przykładzie obciętego normalnego .

Istnieje wiele różnych metod, niektóre proste, a niektóre stosunkowo wydajne.

Zilustruję niektóre podejścia na twoim normalnym przykładzie.

Oto jedna bardzo prosta metoda generowania pojedynczo (w jakimś pseudokodzie):

$\tt{repeat}$ generuj $x_i$ z N (średnia, sd) $\tt{until}$ niższa $\leq x_i\leq$ górna

Jeśli większość rozkładu mieści się w granicach, jest to całkiem rozsądne, ale może być dość powolne, jeśli prawie zawsze generujesz poza granicami.

W R można uniknąć pętli pojedynczej, obliczając obszar w granicach i generując wystarczającą liczbę wartości, dzięki czemu można być niemal pewnym, że po wyrzuceniu wartości poza granice nadal masz tyle wartości, ile potrzeba.
Możesz użyć akceptacji-odrzucenia z jakąś odpowiednią funkcją majorowania w tym przedziale (w niektórych przypadkach mundur będzie wystarczająco dobry). Gdyby granice były dość wąskie w stosunku do SD, ale nie byłeś daleko w tyle, jednolita majoralizacja działałaby dobrze na przykład z normalną.
Jeśli masz dość wydajne cdf i odwrotne cdf (takie jak pnormi qnormdo normalnego rozkładu w R), możesz użyć metody odwrotnej cdf opisanej w pierwszym akapicie sekcji symulacyjnej strony Wikipedii na normalnej skróconej . [W efekcie jest to to samo, co pobranie obciętego munduru (obciętego przy wymaganych kwantylach, który w rzeczywistości nie wymaga żadnych odrzuceń, ponieważ jest to tylko kolejny mundur) i zastosowanie do tego odwrotnego normalnego cdf. Pamiętaj, że może się to nie powieść, jeśli jesteś daleko w ogonie]
Istnieją inne podejścia; ta sama strona Wikipedii wspomina o dostosowaniu metody ziggurat , która powinna działać dla różnych dystrybucji.

Ten sam link do Wikipedii wymienia dwa konkretne pakiety (oba w CRAN) z funkcjami do generowania obciętych normalnych:

MSMOpakowanie R ma funkcję rtnorm, która oblicza czerpie ściętego normalne. truncnormOpakowanie R ma również funkcje czerpania ze ściętego normalne.

Rozglądając się, wiele z tego jest zawartych w odpowiedziach na inne pytania (ale nie do końca duplikowane, ponieważ to pytanie jest bardziej ogólne niż tylko okrojona norma) ... zobacz dodatkową dyskusję w

za. Ta odpowiedź

b. Odpowiedź Xi'ana tutaj , która zawiera link do jego artykułu arXiv (wraz z kilkoma innymi wartościowymi odpowiedziami).

— Glen_b - Przywróć Monikę
źródło

2

Szybkim i brudnym podejściem jest użycie reguły 68-95-99.7 .

W rozkładzie normalnym 99,7% wartości mieści się w 3 standardowych odchyleniach średniej. Tak więc, jeśli ustawisz średnią na żądaną wartość minimalną i maksymalną i ustawisz odchylenie standardowe na 1/3 średniej, otrzymasz (głównie) wartości mieszczące się w pożądanym przedziale. Następnie możesz po prostu posprzątać resztę.

minVal <- 0
maxVal <- 100
mn <- (maxVal - minVal)/2
# Generate numbers (mostly) from min to max
x <- rnorm(count, mean = mn, sd = mn/3)
# Do something about the out-of-bounds generated values
x <- pmax(minVal, x)
x <- pmin(maxVal, x)

Ostatnio napotkałem ten sam problem, próbując wygenerować losowe oceny uczniów na potrzeby danych testowych. W powyższym kodzie użyłem pmaxi pminzamieniłem wartości spoza zakresu na wartość minimalną lub maksymalną. Działa to w moim celu, ponieważ generuję dość małe ilości danych, ale w przypadku większych ilości spowoduje to zauważalne nierówności przy wartościach minimalnej i maksymalnej. Dlatego w zależności od celów lepiej odrzucić te wartości, zastąpić je NAs lub „przerzucić”, dopóki nie znajdą się w granicach.

— Aaron Wells
źródło

Po co to robić? Generowanie zwykłych liczb losowych i upuszczanie tych, które wymagają obcięcia, jest tak proste, że nie trzeba się tym komplikować, chyba że pożądane obcięcie jest bliskie 100% powierzchni gęstości.

— Carl

2

Być może źle interpretuję pierwotne pytanie. Natknąłem się na to pytanie, próbując dowiedzieć się, jak osiągnąć zadanie programowania niezwiązane bezpośrednio ze statystykami w języku R, i dopiero teraz zauważyłem, że ta strona to wymiana stosu statystyk, a nie zmiana stosu programowania. :) W moim przypadku chciałem wygenerować określoną liczbę losowych liczb całkowitych o wartościach od 0 do 100 i chciałem, aby wygenerowane wartości spadały na ładną krzywą dzwonową w tym zakresie. Odkąd to napisałem, zdałem sobie sprawę, że sample(x=min:max, prob=dnorm(...))może to być łatwiejszy sposób.

— Aaron Wells,

@Glen_b Aaron Wells wspomina, sample(x=min:max, prob=dnorm(...))co wydaje się nieco krótsze niż twoja odpowiedź.

— Carl

Pamiętaj jednak, że sample()sztuczka jest przydatna tylko wtedy, gdy próbujesz wybrać losowe liczby całkowite lub inny zestaw dyskretnych, predefiniowanych wartości.

— Aaron Wells

1

$a<b$

$\Phi$ $X_1,...,X_N$ $\mu$ $\sigma^2$ $^\dagger$ $a<b$

X_{i} = μ + σ \cdot Φ^{- 1} (U_{i}) U_{1}, . . ., U_{N} \sim IID U [Φ (\frac{a - μ}{σ}), Φ (\frac{b - μ}{σ})] .

$X_i = \mu + \sigma \cdot \Phi^{-1}(U_i) \quad \quad \quad U_1,...,U_N \sim \text{IID U} \Big[ \Phi \Big( \frac{a-\mu}{\sigma} \Big), \Phi \Big( \frac{b-\mu}{\sigma} \Big) \Big].$

Nie ma wbudowanej funkcji generowania wartości z obciętego rozkładu, ale programowanie tej metody przy użyciu zwykłych funkcji do generowania zmiennych losowych jest banalne. Oto prosta Rfunkcja, rtruncnormktóra implementuje tę metodę w kilku wierszach kodu.

rtruncnorm <- function(N, mean = 0, sd = 1, a = -Inf, b = Inf) {
  if (a > b) stop('Error: Truncation range is empty');
  U <- runif(N, pnorm(a, mean, sd), pnorm(b, mean, sd));
  qnorm(U, mean, sd); }

Jest to wektoryzowana funkcja, która wygeneruje Nlosowe zmienne IID ze skróconego rozkładu normalnego. Łatwo byłoby zaprogramować funkcje dla innych skróconych dystrybucji za pomocą tej samej metody. Nie byłoby też zbyt trudne zaprogramowanie powiązanych funkcji gęstości i kwantylu dla skróconego rozkładu.

$^\dagger$ $\mu$ $\sigma^2$ nie są średnią i wariancją rozkładu obciętego.

— Przywróć Monikę
źródło