Ini adalah generalisasi dari Masalah Ulang Tahun yang terkenal : mengingat orang yang memiliki "ulang tahun" yang acak dan terdistribusi secara seragam di antara sekumpulan kemungkinan , berapa kemungkinan bahwa tidak ada ulang tahun yang dibagikan oleh lebih dari orang?d = 6 m = 20n = 100d= 6m = 20
Perhitungan yang tepat menghasilkan jawaban (untuk menggandakan presisi). Saya akan membuat sketsa teorinya dan memberikan kode untuk umum Waktu asimptotik kode adalah yang membuatnya cocok untuk jumlah ulang tahun yang sangat besar dan memberikan kinerja yang wajar sampai berada dalam ribuan. Pada titik itu, perkiraan Poisson yang dibahas pada Memperluas paradoks ulang tahun menjadi lebih dari 2 orang seharusnya bekerja dengan baik dalam banyak kasus.n , m , d . O ( n 2 log ( d ) ) d n0,267747907805267n , m , d.O ( n2catatan( d) )dn
Penjelasan solusinya
Fungsi menghasilkan probabilitas (pgf) untuk hasil gulungan independen sisi adalahdnd
d- nfn( x1, x2, ... , xd) = d- n( x1+ x2+ ⋯ + xd)n.
Koefisien dalam perluasan multinomial ini memberikan sejumlah cara di mana wajah dapat muncul tepat kali , i e i i = 1 , 2 , … , d .xe11xe22⋯ xeddsayaesayai = 1 , 2 , … , d.
Membatasi minat kami untuk tidak lebih dari penampilan oleh wajah apa pun sama dengan mengevaluasi modulo yang ideal dihasilkan oleh Untuk melakukan evaluasi ini, gunakan Teorema Binomial secara rekursif untuk mendapatkanf n I x m + 1 1 , x m + 1 2 , … , x m + 1 d .mfnsayaxm + 11, xm + 12, ... , xm + 1d.
fn( x1, ... , xd)= ( ( x1+ ⋯ + xr) + ( xr + 1+ xr + 2+ ⋯ + x2 r) )n= ∑k = 0n( nk) (x1+ ⋯ + xr)k( xr + 1+ ⋯ + x2 r)n - k= ∑k = 0n( nk) fk( x1, ... , xr) fn - k( xr + 1, ... , x2 r)
ketika adalah genap. Menulis ( istilah ), kami memilikif ( d ) n = f n ( 1 , 1 , … , 1 ) dd= 2 rf( d)n= fn( 1 , 1 , … , 1 )d
f( 2 r )n= ∑k = 0n( nk) f( r )kf( r )n - k.(Sebuah)
Ketika adalah ganjil, gunakan dekomposisi analogd= 2 r + 1
fn( x1, ... , xd)= ( ( x1+ ⋯ + x2 r) + x2 r + 1)n= ∑k = 0n( nk) fk( x1, ... , x2 r) fn - k( x2 r + 1) ,
memberi
f( 2 r + 1 )n= ∑k = 0n( nk) f( 2 r )kf( 1 )n - k.(b)
Dalam kedua kasus, kita juga dapat mengurangi semuanya modulo , yang mudah dilakukan sejak awalsaya
fn( xj) ≅{ xn0n ≤ mn > mmodsaya,
memberikan nilai awal untuk rekursi,
f( 1 )n= { 10n ≤ mn > m
Apa yang membuat ini efisien adalah bahwa dengan memisahkan variabel menjadi dua kelompok yang berukuran sama dari masing-masing variabel dan mengatur semua nilai variabel menjadi kita hanya perlu mengevaluasi semuanya satu kali untuk satu kelompok dan kemudian menggabungkan hasilnya. Ini membutuhkan komputasi hingga istilah, masing-masing membutuhkan perhitungan untuk kombinasi. Kita bahkan tidak memerlukan array 2D untuk menyimpan , karena ketika menghitung hanya dan yang diperlukan.dr1,n+1O(n)f(r)nf(d)n,f(r)nf(1)n
Jumlah langkah adalah satu kurang dari jumlah digit dalam ekspansi biner dari (yang menghitung pemisahan menjadi kelompok yang sama dalam rumus ) ditambah jumlah yang ada dalam ekspansi (yang menghitung semua kali ganjil) nilai dijumpai, membutuhkan penerapan rumus ). Itu masih hanya langkah .d(a)(b)O(log(d))
Pada R
workstation yang berumur satu dekade, pekerjaan dilakukan dalam 0,007 detik. Kode terdaftar di akhir posting ini. Ia menggunakan logaritma probabilitas, bukan probabilitas itu sendiri, untuk menghindari kemungkinan meluap atau mengakumulasi terlalu banyak aliran. Ini memungkinkan untuk menghapus faktor dalam solusi sehingga kami dapat menghitung jumlah yang mendasari probabilitas.d−n
Perhatikan bahwa prosedur ini menghasilkan komputasi seluruh urutan probabilitas sekaligus, yang dengan mudah memungkinkan kita mempelajari bagaimana peluang berubah dengan .f0,f1,…,fnn
Aplikasi
Distribusi dalam Masalah Ulang Tahun yang digeneralisasi dihitung oleh fungsi tmultinom.full
. Satu-satunya tantangan terletak pada menemukan batas atas untuk jumlah orang yang harus hadir sebelum peluang -collision menjadi terlalu besar. Kode berikut melakukan ini dengan kekerasan, mulai dengan kecil dan menggandakannya sampai cukup besar. Karena itu seluruh perhitungan membutuhkan waktu di mana adalah solusinya. Seluruh distribusi probabilitas untuk jumlah orang melalui dihitung.n O ( n 2 log ( n ) log ( d ) ) n nm+1nO(n2log(n)log(d))nn
#
# The birthday problem: find the number of people where the chance of
# a collision of `m+1` birthdays first exceeds `alpha`.
#
birthday <- function(m=1, d=365, alpha=0.50) {
n <- 8
while((p <- tmultinom.full(n, m, d))[n] > alpha) n <- n * 2
return(p)
}
Sebagai contoh, jumlah minimum orang yang dibutuhkan dalam kerumunan untuk membuatnya lebih mungkin daripada tidak setidaknya delapan dari mereka berbagi ulang tahun adalah , seperti yang ditemukan dalam perhitungan . Hanya perlu beberapa detik. Berikut adalah plot bagian dari output:798birthday(7)
Versi khusus dari masalah ini dibahas di Memperluas paradoks ulang tahun ke lebih dari 2 orang , yang berkenaan dengan kasus dadu sisi yang digulung beberapa kali.365
Kode
# Compute the chance that in `n` independent rolls of a `d`-sided die,
# no side appears more than `m` times.
#
tmultinom <- function(n, m, d, count=FALSE) tmultinom.full(n, m, d, count)[n+1]
#
# Compute the chances that in 0, 1, 2, ..., `n` independent rolls of a
# `d`-sided die, no side appears more than `m` times.
#
tmultinom.full <- function(n, m, d, count=FALSE) {
if (n < 0) return(numeric(0))
one <- rep(1.0, n+1); names(one) <- 0:n
if (d <= 0 || m >= n) return(one)
if(count) log.p <- 0 else log.p <- -log(d)
f <- function(n, m, d) { # The recursive solution
if (d==1) return(one) # Base case
r <- floor(d/2)
x <- double(f(n, m, r), m) # Combine two equal values
if (2*r < d) x <- combine(x, one, m) # Treat odd `d`
return(x)
}
one <- c(log.p*(0:m), rep(-Inf, n-m)) # Reduction modulo x^(m+1)
double <- function(x, m) combine(x, x, m)
combine <- function(x, y, m) { # The Binomial Theorem
z <- sapply(1:length(x), function(n) { # Need all powers 0..n
z <- x[1:n] + lchoose(n-1, 1:n-1) + y[n:1]
z.max <- max(z)
log(sum(exp(z - z.max), na.rm=TRUE)) + z.max
})
return(z)
}
x <- exp(f(n, m, d)); names(x) <- 0:n
return(x)
}
Jawabannya didapat dengan
print(tmultinom(100,20,6), digits=15)
0,267747907805267
Perhitungan brute force
Kode ini membutuhkan beberapa detik di laptop saya
output: 0,2677479
Tetapi tetap mungkin menarik untuk menemukan metode yang lebih langsung jika Anda ingin melakukan banyak perhitungan ini atau menggunakan nilai yang lebih tinggi, atau hanya demi mendapatkan metode yang lebih elegan.
Setidaknya perhitungan ini memberikan angka yang dihitung secara sederhana, tetapi valid, untuk memeriksa metode lain (yang lebih rumit).
sumber