Konversi MCMC ke nilai tunggal?

13

Saya mencoba menyesuaikan model hierarkis menggunakan jags, dan paket rjags. Variabel hasil saya adalah y, yang merupakan urutan uji coba bernoulli. Saya memiliki 38 subjek manusia yang berkinerja di bawah dua kategori: P dan M. Berdasarkan analisis saya, setiap pembicara memiliki probabilitas keberhasilan dalam kategori P dari dan probabilitas keberhasilan dalam kategori M dari . Saya juga berasumsi bahwa ada beberapa hyperparameter tingkat komunitas untuk P dan M: dan .θ p × θ m μ p μ mθpθp×θmμpμm

Jadi, untuk setiap pembicara: dan mana dan mengontrol bagaimana memuncak distribusi sekitar dan .θ mb e t a ( μ m × κ m , ( 1 - μ m ) × κ m ) κ p κ m μ p μ mθpbeta(μp×κp,(1μp)×κp)θmbeta(μm×κm,(1μm)×κm)κpκmμpμm

Juga , .μ mb e t a ( A m , B m )μpbeta(Ap,Bp)μmbeta(Am,Bm)

Inilah model jags saya:

model{
## y = N bernoulli trials
## Each speaker has a theta value for each category
for(i in 1:length(y)){
    y[i] ~ dbern( theta[ speaker[i],category[i]])
}

## Category P has theta Ptheta
## Category M has theta Ptheta * Mtheta
## No observed data for pure Mtheta
##
## Kp and Km represent how similar speakers are to each other 
## for Ptheta and Mtheta
for(j in 1:max(speaker)){
    theta[j,1] ~ dbeta(Pmu*Kp, (1-Pmu)*Kp)
    catM[j] ~ dbeta(Mmu*Km, (1-Mmu)*Km)
    theta[j,2] <- theta[j,1] * catM[j]
}

## Priors for Pmu and Mmu
Pmu ~ dbeta(Ap,Bp)
Mmu ~ dbeta(Am,Bm)

## Priors for Kp and Km
Kp ~ dgamma(1,1/50)
Km ~ dgamma(1,1/50)

## Hyperpriors for Pmu and Mmu
Ap ~ dgamma(1,1/50)
Bp ~ dgamma(1,1/50)
Am ~ dgamma(1,1/50)
Bm ~ dgamma(1,1/50)
}

Masalah yang saya miliki adalah ketika saya menjalankan model ini dengan 5000 iterasi untuk beradaptasi, kemudian mengambil 1000 sampel, Mmudan Kmtelah terkonvergensi menjadi nilai tunggal. Saya telah menjalankannya dengan 4 rantai, dan setiap rantai tidak memiliki nilai yang sama, tetapi di dalam setiap rantai hanya ada satu nilai.

Saya cukup baru untuk menyesuaikan model hierarkis menggunakan metode MCMC, jadi saya bertanya-tanya seberapa buruk ini. Haruskah saya menganggap ini sebagai tanda bahwa model ini tidak sesuai harapan, bahwa ada sesuatu yang salah dengan prior saya, atau apakah ini setara dengan kursus?

Sunting: Dalam kasus itu penting, nilai untuk itu konvergen (rata-rata lintas) adalah 0,91 dan adalah 1,78κ mμmκm

JoFrhwld
sumber
Jika saya memahami Anda dengan benar, parameter ini "menyatu" pada satu nilai tetap di setiap rantai (setelah beberapa iterasi tidak berubah sama sekali), tetapi nilai itu berbeda untuk setiap rantai yang Anda jalankan? Itu terdengar buruk, seperti mungkin langkah Metropolis Hastings yang benar-benar jelek. Bisa jadi model Anda, bisa juga JAGS, bisa jadi kombinasi keduanya. Agaknya model ini tidak terlalu lama untuk dipasangkan, jadi saya akan mencoba menjalankan (lebih banyak) rantai lebih dulu, terutama untuk periode adaptasi.
JMS
Jadi, saya memperbarui model dengan 5000 lebih banyak iterasi, dan parameter yang dipermasalahkan tidak bergerak. Saya tidak menyadari bahwa mereka dapat jatuh ke dalam minimum lokal seperti ini.
JoFrhwld
1
saran cepat: 1. Coba gunakan dbin, dengan n = 1. Dan gunakan batas untuk membatasi nilai p. Sesuatu seperti ini: p.bound [i] <- maks (0, min (1, p [i]))
Manoel Galdino
1
Beberapa pertanyaan klarifikasi: 1. Anda memiliki 38 subjek dalam kategori P dan 38 subyek dalam kategori M, seperti panjang (y) = 76? 2. Bisakah Anda memberikan lebih banyak informasi latar belakang tentang alasan untuk hyperparamters dan eksperimen? Agak membingungkan bagi saya.
Manoel Galdino
1
Saya mungkin hanya memperbaiki semua hyperparameters, seperti dalam theta [j, 1] ~ dbeta (1.1, 1.1) atau sesuatu, dan lihat hasil apa yang Anda dapatkan sebelum mencoba beralih ke hyperprior. Juga, prior theta [j, 2] adalah produk dari dua beta, yang secara umum bukan beta itu sendiri, dan tentu saja theta [j, 2] sebagai hasilnya. Sepertinya Anda yang dimaksudkan ini; hanya karena penasaran, mengapa?
jbowman

Jawaban:

2

Ini lebih merupakan komentar, tetapi karena saya tidak memiliki reputasi yang cukup, saya mungkin juga menjawab.

Dari pengalaman saya yang terbatas dengan MCMC samplers, yang saya amati adalah bahwa parameternya cenderung tetap ketika hyperparameter terlalu sempit. Ketika mereka mengontrol penyebaran parameter, mereka mencegah ruang solusi menjadi sampel yang efisien.

Cobalah untuk mengambil nilai yang lebih besar untuk hyperparameter dan lihat apa yang terjadi.

Makalah teknis ini banyak membantu saya memahami sampler MCMC. Ini terdiri dari dua sampler, Gibbs (yang Anda gunakan), dan Hybrid Monte Carlo, dan dengan cepat menjelaskan cara memilih prior, hyperpriors dan nilai untuk parameter dan hyperparameters.

Giezi
sumber
0

Ini bisa menjadi masalah struktur rantai. Di mana Anda berakhir tergantung di mana Anda mulai. Untuk menggunakan MCMC Anda ingin rantai itu berulang yang berarti bahwa di mana pun Anda memulai, Anda bisa mencapai setiap negara bagian di ruang keadaan. Jika rantai tidak berulang Anda dapat terjebak dalam subset ruang keadaan. Gagasan MCMC adalah untuk memiliki distribusi stasioner yang sudah ada yang pada akhirnya rantai akan berakhir. Distribusi stasioner ini biasanya memiliki probabilitas positif untuk berada di negara bagian mana pun dalam rantai dan tidak berakhir terperangkap pada satu titik seperti yang Anda jelaskan . Saya tidak dapat memeriksa algoritme Anda, tetapi mungkin Anda memiliki kesalahan di dalamnya. Mungkin juga Anda telah mendefinisikan masalah di mana rantai Markov Anda tidak berulang.

Jika Anda ingin berpengetahuan luas di MCMC, saya sarankan Anda membaca Buku Pegangan Markov Chain Monte Carlo yang memiliki artikel yang menggambarkan setiap aspek MCMC. Itu diterbitkan oleh CRC Press pada tahun 2011.

Michael R. Chernick
sumber