Jak próbkować z dyskretnego rozkładu na liczbach całkowitych nieujemnych?

10

Mam następujący dyskretny rozkład, w którym są znanymi stałymi:α,β

p(x;α,β)=Beta(α+1,β+x)Beta(α,β)dla x=0,1,2),

Jakie są podejścia do skutecznego próbkowania z tej dystrybucji?

jII
źródło

Odpowiedzi:

9

Jest to rozkład dwumianowy Beta-ujemny z parametrem w twoim przypadku, z wykorzystaniem notacji Wikipedia. Nazwano także rozkładem Beta-Pascal, gdy jest liczbą całkowitą. Jak zauważyłeś w komentarzu, jest to rozkład predykcyjny w bayesowskim ujemnym modelu dwumianowym ze sprzężoną Beta przed prawdopodobieństwem sukcesu.r=1r

W ten sposób możesz próbkować go, próbkując zmienną a następnie próbkując ujemną zmienną dwumianową (przy twoim przypadku , to znaczy powiedzieć rozkład geometryczny).Beta(α,β)uNB(r,u)r=1

Rozkład ten jest realizowany w pakiecie R brr. Próbnik ma nazwę rbeta_nbinom, pmf ma nazwę dbeta_nbinomitp. Notacje to , , . Czek:a=rc=αd=β

> Alpha <- 2; Beta <- 3
> a <- 1
> all.equal(brr::dbeta_nbinom(0:10, a, Alpha, Beta), beta(Alpha+a, Beta+0:10)/beta(Alpha,Beta))
[1] TRUE

Patrząc na kod, widać, że faktycznie wywołuje on ghyper(uogólnioną hipergeometryczną) rodzinę dystrybucji SuppDistspakietu:

brr::rbeta_nbinom
function(n, a, c, d){
  rghyper(n, -d, -a, c-1)
}

Innymi słowy, rozkład BNB jest znany jako uogólniony rozkład hipergeometryczny typu IV . Zobacz pomoc ghyperw SuppDistspakiecie. Wierzę, że można to również znaleźć w książce Johnson & al Univariate Discrete Distribution .

Stéphane Laurent
źródło
Ta odpowiedź jest świetna, ale byłoby jeszcze lepiej, gdybyś udowodnił, że zaksięgowana gęstość OP jest taka sama jak ujemna gęstość dwumianowa.
Sycorax mówi Przywróć Monikę
1
@ user777 Myślę, że autor OP sam to udowodnił, biorąc pod uwagę jego komentarz do odpowiedzi Xiana (tylna predykcyjna dystrybucja w ujemnym modelu dwumianowym z koniugatem Beta przed).
Stéphane Laurent,
10

Jeśli się uwzględni

Beta(α+1,β+x)Beta(α,β)=αα+β+xβ+x1α+β+x1βα+β
zmniejsza się z x, Proponuję wygenerować zmienną jednolitą uU(0,1) i obliczanie skumulowanych kwot
S.k=x=0kBeta(α+1,β+x)Beta(α,β)
aż do
S.k>u
Realizacja jest wtedy równa odpowiadającej k. Od
Rx=Beta(α+1,β+x)Beta(α,β)=αα+β+xβ+x-1α+β+x-1βα+β=α+β+x-1α+β+xβ+x-1α+β+x-1Rx-1=β+x-1α+β+xRx-1
i
S.k=S.k-1+Rk
w obliczeniach można całkowicie uniknąć korzystania z funkcji gamma.
Xi'an
źródło
1
(+1) Korzystanie S.k=1-Γ(za+b)Γ(b+k+1)Γ(b)Γ(za+b+k+1)znacznie przyspieszy pracę.
whuber
1
Re edycja: Podejrzewam, że wykorzystanie funkcji gamma będzie jednak pomocne w rozwiązywaniu problemu k pod względem u, α, i β. Na przykład można znaleźć wstępne przybliżenie dou za pomocą wzoru Stirlinga w ocenie Γ(b+k+1) i Γ(za+b+k+1)a następnie dopracowując to za pomocą kilku kroków Newtona-Raphsona. Potrzebują one oceny log Gamma i jej pochodnej. Oczywiście jeśliα i βoba są integralne, to rozwiązanie jest źródłem wielomianu - ale nawet wtedy używanie Gammy może być nadal dobrym rozwiązaniem.
whuber
1
Świetna odpowiedź! Zaakceptowałem odpowiedź udzieloną przez SL, ponieważ zwróciło mi ona uwagę na kluczową kwestię (nie jest częścią pierwotnego pytania), że próbkowanie z predykcji tylnej jest równoważne próbkowaniu parametru z tylnej, a następnie próbkowanie danych z prawdopodobieństwa. W szczególności powyższa funkcja rozkładu jest tylną predykcją danych geometrycznych z Beta przed parametremp.
jII