Bagaimana saya (secara numerik) dapat memperkirakan nilai untuk distribusi beta dengan alfa & beta besar

Apakah ada cara yang stabil secara numerik untuk menghitung nilai distribusi beta untuk alpha integer besar, beta (misal alpha, beta> 1000000)?

Sebenarnya, saya hanya perlu interval kepercayaan 99% di sekitar mode, jika itu entah bagaimana membuat masalah lebih mudah.

Tambahkan : Maaf, pertanyaan saya tidak sejelas yang saya kira. Yang ingin saya lakukan adalah ini: Saya memiliki mesin yang memeriksa produk pada ban berjalan. Sebagian kecil dari produk ini ditolak oleh mesin. Sekarang jika operator mesin mengubah beberapa pengaturan inspeksi, saya ingin menunjukkan kepadanya perkiraan laju penolakan dan beberapa petunjuk tentang seberapa andal perkiraan tersebut.

Jadi saya pikir saya memperlakukan tingkat penolakan aktual sebagai variabel acak X, dan menghitung distribusi probabilitas untuk variabel acak berdasarkan jumlah objek yang ditolak N dan objek yang diterima M. Jika saya mengasumsikan distribusi seragam sebelumnya untuk X, ini adalah distribusi beta tergantung pada N dan M. Saya bisa menampilkan distribusi ini kepada pengguna secara langsung atau menemukan interval [l, r] sehingga tingkat tolak aktual dalam interval ini dengan p> = 0,99 (menggunakan terminologi shabbychef) dan menampilkan ini selang. Untuk M kecil, N (yaitu segera setelah perubahan parameter), saya dapat menghitung distribusi secara langsung dan mendekati interval [l, r]. Tetapi untuk M besar, N, pendekatan naif ini mengarah pada kesalahan underflow, karena x ^ N * (1-x) ^ M adalah kecil untuk diwakili sebagai pelampung presisi ganda.

Saya kira taruhan terbaik saya adalah menggunakan beta-distribusi naif saya untuk M kecil, N dan beralih ke distribusi normal dengan rata-rata dan varians yang sama begitu M, N melebihi ambang batas tertentu. Apakah itu masuk akal?

Perkiraan normal bekerja dengan sangat baik, terutama di bagian ekor. Gunakan rata-rata dan varian α $\alpha/(\alpha+\beta)$ . Misalnya, kesalahan relatif absolut dalam probabilitas ekor dalam situasi sulit (di mana kemiringan mungkin menjadi perhatian) sepertiα=106,β=108puncak di sekitar0,00026dan kurang dari0,00006ketika Anda lebih dari 1 SD dari mean. (Inibukankarena beta sangat besar: denganα=β=106, kesalahan relatif absolut dibatasi oleh$\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$.) Dengan demikian, perkiraan ini sangat baik untuk tujuan apa pun yang melibatkan interval 99%.

Sehubungan dengan pengeditan pada pertanyaan, perhatikan bahwa seseorang tidak menghitung beta integral dengan benar-benar mengintegrasikan integrand: tentu saja Anda akan mendapatkan arus bawah (meskipun mereka tidak terlalu penting, karena mereka tidak berkontribusi cukup besar pada integral) . Ada banyak, banyak cara untuk menghitung integral atau perkiraannya, sebagaimana didokumentasikan dalam Johnson & Kotz (Distribusi dalam Statistik). Kalkulator daring ditemukan di http://www.danielsoper.com/statcalc/calc37.aspx . Anda sebenarnya membutuhkan kebalikan dari integral ini. Beberapa metode untuk menghitung kebalikannya didokumentasikan di situs Mathematica di http://functions.wolfram.com/GammaBetaErf/InverseBetaRegularized/. Kode disediakan dalam Numerical Recipes (www.nr.com). Kalkulator online yang sangat bagus adalah situs Wolfram Alpha (www.wolframalpha.com): masukkan inverse beta regularized (.005, 1000000, 1000001)untuk titik akhir kiri dan inverse beta regularized (.995, 1000000, 1000001)untuk titik akhir kanan ( , interval 99%).$\alpha=1000000, \beta=1000001$

Eksperimen grafis cepat menunjukkan bahwa distribusi beta terlihat sangat seperti distribusi normal ketika alfa dan beta keduanya sangat besar. Dengan googling "batas distribusi beta normal" saya menemukan http://nrich.maths.org/discus/messages/117730/143065.html?1200700623 , yang memberikan 'bukti' handwaving.

Halaman wikipedia untuk distribusi beta memberikan rata-rata, mode (v dekat ke rata-rata untuk alpha dan beta besar) dan varians, sehingga Anda dapat menggunakan distribusi normal dengan rata-rata & varians yang sama untuk mendapatkan perkiraan. Apakah itu perkiraan yang cukup baik untuk tujuan Anda tergantung pada apa tujuan Anda.

$[l,r]$$l$$r$$[l,r]$$\alpha, \beta$ $l$$r$ sebagai angka yang berbeda, jadi rute ini mungkin cukup baik.

$p$$p$$log(p/(1-p))$$min(\alpha,\beta) > 100$

Sebagai contoh

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)))

biasanya menghasilkan keluaran seperti

ringkasan (replikasi (50, f (10000, 100, 1000000))) Min. 1 Qu. Median Mean 3 Qu. Maks. 0,01205 0,10870 0,18680 0,24810 0,36170 0,68730

yaitu nilai-p khas sekitar 0,2.

$\alpha=100, \beta=100000$

$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)))

menghasilkan sesuatu seperti

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 

dengan nilai-p khas sekitar 0,01

Fungsi R qqnormjuga memberikan visualisasi yang bermanfaat, menghasilkan plot yang sangat lurus untuk distribusi log-odds yang menunjukkan perkiraan normalitas distribusi variabel beta dsitribute menghasilkan kurva khas yang menunjukkan non normalitas.

$\alpha,\beta$