Jak mogę (liczbowo) aproksymować wartości dla rozkładu beta z dużymi wartościami alfa i beta

Czy istnieje stabilny numerycznie sposób obliczania wartości rozkładu beta dla dużej liczby całkowitej alfa, beta (np. Alfa, beta> 1000000)?

Właściwie potrzebuję tylko 99% przedziału ufności wokół trybu, jeśli to w jakiś sposób ułatwi problem.

Dodaj : Przepraszam, moje pytanie nie było tak jasno określone, jak myślałem. Chcę to zrobić: mam maszynę, która sprawdza produkty na przenośniku taśmowym. Część frakcji tych produktów jest odrzucana przez maszynę. Teraz, jeśli operator maszyny zmieni jakieś ustawienia inspekcji, chcę mu pokazać szacunkową częstotliwość odrzucania i podpowiedź na temat tego, jak wiarygodne jest bieżące oszacowanie.

Pomyślałem więc, że traktuję rzeczywistą częstotliwość odrzucania jako zmienną losową X i obliczam rozkład prawdopodobieństwa dla tej zmiennej losowej na podstawie liczby odrzuconych obiektów N i zaakceptowanych obiektów M. Jeśli założę jednolity wcześniejszy rozkład dla X, jest to rozkład beta w zależności od N i M. Mogę albo wyświetlić ten rozkład użytkownikowi bezpośrednio, albo znaleźć przedział [l, r], aby rzeczywista częstość odrzucania była w tym przedziale przy p> = 0,99 (używając terminologii shabbychef) i wyświetlić to interwał. Dla małych M, N (tj. Bezpośrednio po zmianie parametru) mogę obliczyć rozkład bezpośrednio i przybliżać przedział [l, r]. Ale w przypadku dużych M, N takie naiwne podejście prowadzi do błędów niedopełnienia, ponieważ x ^ N * (1-x) ^ M jest zbyt małe, aby można je było przedstawić jako zmiennoprzecinkowe podwójnej precyzji.

Myślę, że moim najlepszym wyborem jest użycie mojej naiwnej dystrybucji beta dla małych M, N i przejście do normalnej dystrybucji z tą samą średnią i wariancją, gdy tylko M, N przekroczy pewien próg. Czy to ma sens?

Normalne przybliżenie działa wyjątkowo dobrze, szczególnie w ogonach. Użyj średniej z i wariantu . Na przykład bezwzględny błąd względny w prawdopodobieństwie ogona w trudnej sytuacji (gdzie skośność może budzić obawy), taki jak osiąga około i jest mniejszy niż gdy jesteś więcej niż 1 SD od średniej. ( Nie dzieje się tak, ponieważ beta jest tak duża: przy bezwzględne błędy względne są ograniczone przez$\alpha/(\alpha+\beta)$$\frac{\alpha\beta}{(\alpha+\beta)^{2} (1+\alpha+\beta)}$$\alpha = 10^6, \beta = 10^8$$0.00026$$0.00006$$\alpha = \beta = 10^6$$0.0000001$.) To przybliżenie jest zatem doskonałe do zasadniczo dowolnego celu obejmującego przedziały 99%.

W świetle zmian wprowadzonych do pytania należy zauważyć, że nie oblicza się całek beta poprzez faktyczne całkowanie całki: oczywiście dostaniesz niedopełnienia (chociaż tak naprawdę nie mają znaczenia, ponieważ nie mają znaczącego wpływu na całkę) . Istnieje wiele sposobów obliczania całki lub przybliżania jej, jak udokumentowano w Johnson & Kotz (Distribution in Statistics). Kalkulator online znajduje się na stronie http://www.danielsoper.com/statcalc/calc37.aspx . W rzeczywistości potrzebujesz odwrotności tej całki. Niektóre metody obliczania odwrotności są udokumentowane na stronie Mathematica pod adresem http://functions.wolfram.com/GammaBetaErf/InverseBetaRegularized/. Kod znajduje się w przepisach numerycznych (www.nr.com). Naprawdę fajnym kalkulatorem online jest strona Wolfram Alpha (www.wolframalpha.com): wprowadź inverse beta regularized (.005, 1000000, 1000001)lewy punkt końcowy i inverse beta regularized (.995, 1000000, 1000001)prawy punkt końcowy ( , przedział 99%).$\alpha=1000000, \beta=1000001$

Szybki eksperyment graficzny sugeruje, że rozkład beta wygląda bardzo podobnie do rozkładu normalnego, gdy alfa i beta są bardzo duże. Po przejrzeniu „normalnego limitu dystrybucji beta” znalazłem http://nrich.maths.org/discus/messages/117730/143065.html?1200700623 , co daje „dowód” na machanie ręką.

Strona wikipedii dla dystrybucji beta podaje jej średnią, tryb (v zbliżony do średniej dla dużej alfa i beta) i wariancję, więc można użyć rozkładu normalnego z tą samą średnią i wariancją, aby uzyskać przybliżenie. To, czy jest to wystarczająco dobre przybliżenie dla twoich celów, zależy od twoich celów.

Mam zamiar wywnioskować, że chcesz przedział taki, że prawdopodobieństwo losowego losowania z Beta RV jest w przedziale z prawdopodobieństwem 0,99, z dodatkowymi punktami dla i symetrycznymi wokół trybu. Na podstawie Nierówności Gaussa lub nierówności Vysochanskii-Petunin można konstruować przedziały, które zawierają przedział , i byłyby to całkiem przyzwoite przybliżenia. W przypadku dostatecznie dużych będziesz miał problemy z niedopełnieniem liczb, nawet reprezentując i jako odrębne liczby, więc ta droga może być wystarczająco dobra.$[l,r]$$l$$r$$[l,r]$$\alpha, \beta$ $l$$r$

Jeśli jest zmienną rozproszoną w wersji beta, to logarytmiczna szansa (tj. jest w przybliżeniu normalnie rozłożona. Jest to prawdą nawet w przypadku mocno wypaczonych rozkładów beta wraz zp l o g ( p / ( 1 - p ) ) m i n ( α , β ) > $p$$p$$log(p/(1-p))$$min(\alpha,\beta) > 100$

Na przykład

f <- function(n, a, b) {
    p <- rbeta(n, a, b)
    lor <- log(p/(1-p))
    ks.test(lor, 'pnorm', mean(lor), sd(lor))$p.value
}
summary(replicate(50, f(10000, 100, 1000000)))

zazwyczaj daje wynik podobny do

podsumowanie (powtórzenie (50, f (10000, 100, 1000000))) min. 1st Qu. Mediana oznacza 3. kwartę Max. 0,01205 0,108 0,08680 0,248 0,36670 0,68730

tzn. typowe wartości p wynoszą około 0,2.

Tak więc nawet przy 10000 próbkach test Kołmogorowa-Smirnowa nie ma mocy, aby odróżnić transformację ilorazu logarytmicznego silnie wypaczonej zmiennej rozproszonej beta o .$\alpha=100, \beta=100000$

Jednak podobny test rozkładu samego$p$

f2 <- function(n, a, b) {
    p <- rbeta(n, a, b)
    ks.test(p, 'pnorm', mean(p), sd(p))$p.value
}
summary(replicate(50, f2(10000, 100, 1000000)))

produkuje coś podobnego

summary(replicate(50, f2(10000, 100, 1000000)))
     Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
2.462e-05 3.156e-03 7.614e-03 1.780e-02 1.699e-02 2.280e-01 

z typowymi wartościami p około 0,01

Funkcja R qqnormdaje również pomocną wizualizację, tworząc bardzo prosto wyglądający wykres rozkładu logarytmiczno-szansowego wskazujący przybliżoną normalność, rozkład zmiennej dsitribute beta daje charakterystyczną krzywą wskazującą na nienormalność

Dlatego rozsądne jest zastosowanie aproksymacji Gaussa w przestrzeni logarytmicznej, nawet dla mocno wypaczonych wartości o ile oba są większe niż 100.$\alpha,\beta$