Pemecahan sampel secara analitik dengan atau tanpa penggantian setelah Poisson / Binomial negatif

Versi pendek

Saya mencoba untuk secara analitis menyelesaikan / memperkirakan kemungkinan gabungan yang dihasilkan dari Poisson independen menarik dan pengambilan sampel lebih lanjut dengan atau tanpa penggantian (saya tidak benar-benar peduli yang mana). Saya ingin menggunakan kemungkinan dengan MCMC (Stan), jadi saya perlu solusinya hanya sampai jangka waktu yang konstan. Pada akhirnya, saya ingin membuat model proses di mana penarikan awal dari neg. distribusi binomial, tapi saya pikir saya akan bisa sampai di sana dengan solusi untuk kasus Poisson.

Sangat mungkin bahwa solusinya tidak layak (saya tidak mengerti matematika cukup untuk dapat mengetahui apakah ini masalah yang sederhana atau sangat sulit). Karena itu saya juga tertarik pada perkiraan, hasil negatif atau intuisi mengapa masalahnya mungkin tidak dapat diselesaikan (misalnya membandingkan dengan masalah sulit yang diketahui). Tautan ke makalah / teorema / trik yang berguna yang akan membantu saya bergerak maju adalah jawaban yang baik meskipun koneksi mereka ke masalah yang ada belum sepenuhnya teratasi.

Pernyataan formal

Lebih formal, pertama diambil secara independen dan kemudian saya sampel item secara acak dari semua untuk mendapatkan . Yaitu saya menggambar bola berwarna dari guci di mana jumlah bola warna diambil dari . Di sini, diasumsikan diketahui dan diperbaiki dan kami mengkondisikan pada . Secara teknis pengambilan sampel dilakukan tanpa penggantian, tetapi mengasumsikan pengambilan sampel dengan penggantian tidak boleh menjadi masalah besar.$Y = (y_1, ..., y_N), y_n \sim Pois(\lambda_n)$$k$$Y$$Z = (z_1,...,z_N)$$k$$n$$Pois(\lambda_n)$$k$$\sum_n y_n \geq k$

Saya telah mencoba dua pendekatan untuk menyelesaikan pengambilan sampel tanpa penggantian (karena ini tampaknya menjadi kasus yang lebih mudah karena beberapa ketentuan membatalkan), tetapi terjebak dengan keduanya. Kemungkinan pengambilan sampel tanpa penggantian adalah:

EDIT: Bagian "percobaan solusi telah dihapus karena solusi dalam jawaban tidak dibangun di atasnya (dan jauh lebih baik)

Kasing tanpa penggantian

Jika Anda memiliki variabel terdistribusi Poisson independen$n$

dan kondisi pada

kemudian

Jadi Anda bisa mengisi guci Anda dengan kali bola berwarna seperti pertama menggambar nilai total (yang merupakan cutoff terdistribusi Poisson oleh kondisi ) dan kemudian mengisi guci dengan bola menurut distribusi multinomial.$n$$Y_i$$K$$K \geq k$$K$

Pengisian guci ini dengan bola , menurut distribusi multinomial, setara dengan menggambar untuk setiap bola secara independen warna dari distribusi kategorikal. Kemudian Anda dapat mempertimbangkan bola pertama yang telah ditambahkan ke guci sebagai mendefinisikan sampel acak (ketika sampel ini diambil tanpa penggantian) dan distribusi untuk ini hanyalah vektor terdistribusi multinomial lain: $K$$k$$\lbrace Z_i \rbrace$

simulasi

##### simulating sample process for 3 variables #######


# settings
set.seed(1)
k = 10
lambda = c(4, 4, 4)
trials = 1000000

# observed counts table
Ocounts = array(0, dim = c(k+1, k+1, k+1))

for (i in c(1:trials)) {
  # draw poisson with limit sum(n) >= k
  repeat {
    Y = rpois(3,lambda)
    if (sum(Y) >= k) {break}
  }
  # setup urn
  urn <- c(rep(1, Y[1]), rep(2, Y[2]), rep(3, Y[3]))
  # draw from urn
  s <- sample(urn, k, replace=0)
  Z = c(sum(s==1),sum(s==2),sum(s==3))
  Ocounts[Z[1]+1, Z[2]+1, Z[3]+1] = Ocounts[Z[1]+1, Z[2]+1, Z[3]+1] + 1
}



# comparison
observed = rep(0, 0.5*k*(k+1))
expected = rep(0, 0.5*k*(k+1))   
t = 1

for (z1 in c(0:k)) {
  for (z2 in c(0:(k-z1))) {  
    z3 = k-z1-z2 
    observed[t] = Ocounts[z1+1, z2+1, z3+1]
    expected[t] = trials*dmultinom(c(z1, z2, z3), prob=lambda)
    t = t+1
  }
}

plot(observed,expected)
x2 <- sum((observed-expected)^2/expected)
pvalue <- 1-pchisq(x2, 66-1)

hasil

> # results from chi-sq test
> x2
[1] 75.49286
> pvalue
[1] 0.1754805

Binomial negatif

Argumen akan bekerja sama untuk kasus distribusi binomial negatif yang menghasilkan ( dalam kondisi tertentu ) menjadi distribusi Dirichlet-multinomial.

Di bawah ini adalah contoh simulasi

##### simulating sample process for 3 variables #######

# dirichlet multinomial for vectors of size 3
ddirmultinom =  function(x1,x2,x3,p1,p2,p3) {
  (factorial(x1+x2+x3)*gamma(p1+p2+p3)/gamma(x1+x2+x3+p1+p2+p3))/
  (factorial(x1)*gamma(p1)/gamma(x1+p1))/
  (factorial(x2)*gamma(p2)/gamma(x2+p2))/
  (factorial(x3)*gamma(p3)/gamma(x3+p3))
}

# settings
set.seed(1)
k = 10
theta = 1
lambda = c(4,4,4)
trials = 1000000

# calculating negative binomials pars
means = lambda
vars = lambda*(1+theta)

ps = (vars-means)/(vars)
rs = means^2/(vars-means)


# observed counts table
Ocounts = array(0, dim = c(k+1,k+1,k+1))

for (i in c(1:trials)) {
  # draw poisson with limit sum(n) >= k
  repeat {
    Y = rnbinom(3,rs,ps)
    if (sum(Y) >= k) {break}
  }
  # setup urn
  urn <- c(rep(1,Y[1]),rep(2,Y[2]),rep(3,Y[3]))
  # draw from urn
  s <- sample(urn,k,replace=0)
  Z = c(sum(s==1),sum(s==2),sum(s==3))
  Ocounts[Z[1]+1,Z[2]+1,Z[3]+1] = Ocounts[Z[1]+1,Z[2]+1,Z[3]+1] + 1
}



# comparison
observed = rep(0,0.5*k*(k+1))
expected = rep(0,0.5*k*(k+1))   
t = 1

for (z1 in c(0:k)) {
  for (z2 in c(0:(k-z1))) {  
    z3 = k-z1-z2 
    observed[t]=Ocounts[z1+1,z2+1,z3+1]
    expected[t]=trials*ddirmultinom(z1,z2,z3,lambda[1]/theta,lambda[2]/theta,lambda[3]/theta)
    t = t+1
  }
}

plot(observed,expected)
x2 <- sum((observed-expected)^2/expected)
pvalue <- 1-pchisq(x2,66-1)

# results from chi-sq test
x2
pvalue