Bagaimana kabar Garvan?
Masalahnya adalah kita tidak tahu berapa jumlah nol yang diamati. Kami harus memperkirakan ini. Prosedur statistik klasik untuk situasi seperti ini adalah algoritma Ekspektasi-Maksimalisasi.
Contoh sederhana:
Asumsikan kita menarik dari populasi yang tidak diketahui (1.000.000) dengan konstanta poisson 0,2.
counts <- rpois(1000000, 0.2)
table(counts)
0 1 2 3 4 5
818501 164042 16281 1111 62 3
Tapi kami tidak mengamati jumlah nol. Sebaliknya, kami mengamati ini:
table <- c("0"=0, table(counts)[2:6])
table
0 1 2 3 4 5
0 164042 16281 1111 62 3
Frekuensi yang mungkin diamati
k <- c("0"=0, "1"=1, "2"=2, "3"=3, "4"=4, "5"=5)
Inisialisasi mean distribusi Poisson - cukup tebak (kami tahu 0,2 di sini).
lambda <- 1
Harapan - Distribusi Poisson
P_k <- lambda^k*exp(-lambda)/factorial(k)
P_k
0 1 2 3 4 5
0.367879441 0.367879441 0.183939721 0.061313240 0.015328310 0.003065662
n0 <- sum(table[2:6])/(1 - P_k[1]) - sum(table[2:6])
n0
0
105628.2
table[1] <- 105628.2
Maksimalisasi
lambda_MLE <- (1/sum(table))*(sum(table*k))
lambda_MLE
[1] 0.697252
lambda <- lambda_MLE
Iterasi kedua
P_k <- lambda^k*exp(-lambda)/factorial(k)
n0 <- sum(table[2:6])/(1 - P_k[1]) - sum(table[2:6])
table[1] <- n0
lambda <- (1/sum(table))*(sum(table*k))
population lambda_MLE
[1,] 361517.1 0.5537774
Sekarang beralih hingga konvergensi:
for (i in 1:200) {
P_k <- lambda^k*exp(-lambda)/factorial(k)
n0 <- sum(table[2:6])/(1 - P_k[1]) - sum(table[2:6])
table[1] <- n0
lambda <- (1/sum(table))*(sum(table*k))
}
cbind( population = sum(table), lambda_MLE)
population lambda_MLE
[1,] 1003774 0.1994473
Perkiraan populasi kami adalah 1003774 dan laju poisson kami diperkirakan sebesar 0,1994473 - ini adalah proporsi estimasi populasi sampel. Masalah utama yang akan Anda miliki dalam masalah biologis tipikal yang Anda hadapi adalah asumsi bahwa tingkat poisson adalah konstan.
Maaf untuk posting yang bertele-tele - wiki ini tidak benar-benar cocok untuk kode R.
Ini terdengar seperti bentuk 'tanda dan tangkap kembali' alias 'tangkap-tangkap kembali', teknik yang terkenal di bidang ekologi (dan beberapa bidang lain seperti epidemiologi). Bukan wilayah saya, tetapi artikel Wikipedia tentang tanda dan penangkapan kembali terlihat masuk akal, meskipun situasi Anda bukanlah situasi yang dijelaskan oleh metode Lincoln – Petersen di sana.
Saya pikir shabbychef adalah salah satu jalur yang tepat untuk situasi Anda, tetapi menggunakan distribusi Poisson untuk memperkirakan binomial mungkin akan membuat segalanya lebih sederhana dan harus menjadi perkiraan yang sangat baik jika ukuran populasi sangat besar, seperti pada contoh Anda. Saya pikir mendapatkan ekspresi eksplisit untuk perkiraan kemungkinan maksimum ukuran populasi harus cukup mudah (lihat misalnya Wikipedia lagi ), meskipun saya tidak punya waktu untuk menyelesaikan rinciannya sekarang.
sumber
sumber