Apakah metode pengambilan sampel "pentingnya Gibbs" berfungsi?

Saya curiga ini adalah pertanyaan yang tidak biasa dan bersifat eksploratif, jadi tolong ajukan pertanyaan ini kepada saya.

Saya bertanya-tanya apakah seseorang dapat menerapkan gagasan pentingnya pengambilan sampel untuk pengambilan sampel Gibbs. Inilah yang saya maksud: dalam pengambilan sampel Gibbs, kami mengubah nilai satu variabel (atau blok variabel) pada suatu waktu, pengambilan sampel dari probabilitas bersyarat mengingat variabel yang tersisa.

Namun, mungkin tidak mungkin atau tidak mudah untuk mengambil sampel dari probabilitas kondisional yang tepat. Jadi alih-alih, kami mengambil sampel dari distribusi proposal dan menggunakan, misalnya, Metropolis-Hastings (MH).$q$

Sejauh ini bagus. Tapi di sini ada jalan yang berbeda: apa yang terjadi jika kita, alih-alih menggunakan MH, menggunakan gagasan yang sama dengan yang digunakan dalam sampel kepentingan, yaitu kita mengambil sampel dari dan mempertahankan bobot penting dari sampel saat ini?$q$$p/q$

Secara lebih rinci: anggap kita memiliki variabel dan distribusi dengan faktor sehingga . Kami menjaga probabilitas proposal digunakan untuk sampel nilai saat ini dari setiap variabel . Pada setiap langkah kami mengubah subset dari variabel dan memperbarui (hanya faktor dan yang terpengaruh). Kami mengambil sampel dan bobot pentingnya untuk menghitung statistik apa pun yang kami minati.$x_1,\dots,x_n$$\phi_1,\dots,\phi_m$$p \propto \prod_{i=1}^m \phi_i$$q_i$$x_i$$p(x)/q(x)$$p$$q$

Apakah algoritma ini benar? Jika tidak, ada alasan jelas mengapa tidak? Secara intuitif masuk akal bagi saya karena tampaknya melakukan hal yang sama dengan kepentingan sampel, tetapi dengan sampel yang tergantung.

Saya menerapkan ini untuk model jalan acak Gaussian dan mengamati bahwa bobot menjadi lebih kecil dan lebih kecil (tetapi tidak monoton), sehingga sampel awal pada akhirnya memiliki terlalu banyak kepentingan dan mendominasi statistik. Saya cukup yakin implementasinya tidak bermasalah, karena pada setiap langkah saya membandingkan bobot yang diperbarui dengan perhitungan brute force secara eksplisit. Perhatikan bahwa bobot tidak turun tanpa batas ke nol, karena mereka mana dan adalah produk dengan jumlah kepadatan terbatas, dan setiap sampel diperoleh dari distribusi Normal yang jarang jarang menjadi nol.$p/q$$p$$q$

Jadi saya mencoba memahami mengapa bobot turun seperti itu, dan apakah ini merupakan konsekuensi dari metode ini sebenarnya tidak benar.

---

Berikut definisi algoritma yang lebih tepat, yang diterapkan pada jalan acak Gaussian pada variabel . Kode berikut di bawah ini.$X_1,\dots,X_n$

Modelnya hanyalah , dengan tetap pada .$X_i \sim \mathcal N(X_{i-1}, \sigma^2), i = 1,\dots,n$$X_0$$0$

Berat sampel saat ini adalah , di mana adalah kepadatan Gaussian dan adalah distribusi dari mana nilai saat ini telah dijadikan sampel. Awalnya, kami hanya sampel nilai-nilai secara maju, jadi dan bobot awal adalah .$\frac{\prod_i p(x_i)}{\prod_i q(x_i)}$$p$$q$$q = p$$1$

Kemudian pada setiap langkah saya memilih untuk diubah. Saya mencicipi nilai baru untuk dari , jadi kepadatan ini menjadi distribusi proposal yang digunakan baru untuk .$j \in \{1,\dots,n\}$$x'_j$$X_j$$\mathcal N(X_{j-1},\sigma^2)$$X_j$

Untuk memperbarui bobot, saya membaginya dengan kepadatan dan dari nilai lama menurut dan , dan kalikan dengan densitas dan dari nilai baru menurut dan . Ini memperbarui numerator dari bobot.$p(x_j | x_{j-1})$$p(x_{j+1} | x_j)$$x_j$$x_{j-1}$$x_{j+1}$$p(x'_j | x_{j-1})$$p(x_{j+1} | x'_j)$$x'_j$$x_{j-1}$$x_{j+1}$$p$

Untuk memperbarui penyebut , saya gandakan bobotnya dengan proposal lama (dengan demikian menghapusnya dari penyebut) dan membaginya dengan .$q$$q(x_j)$$q(x'_j)$

(Karena saya mencicipi dari normal yang berpusat pada , selalu sama dengan sehingga mereka membatalkan dan implementasinya tidak tidak benar-benar menggunakannya).$x'_j$$x_{j-1}$$q(x'_j)$$p(x'_j | x_{j-1})$

Seperti yang saya sebutkan sebelumnya, dalam kode saya membandingkan perhitungan bobot tambahan ini dengan perhitungan eksplisit yang sebenarnya hanya untuk memastikan.

---

Ini kode untuk referensi.

println("Original sample: " + currentSample);
int flippedVariablesIndex = 1 + getRandom().nextInt(getVariables().size() - 1);
println("Flipping: " + flippedVariablesIndex);
double oldValue = getValue(currentSample, flippedVariablesIndex);
NormalDistribution normalFromBack = getNormalDistribution(getValue(currentSample, flippedVariablesIndex - 1));
double previousP = normalFromBack.density(oldValue);
double newValue = normalFromBack.sample();
currentSample.set(getVariable(flippedVariablesIndex), newValue);
double previousQ = fromVariableToQ.get(getVariable(flippedVariablesIndex));
fromVariableToQ.put(getVariable(flippedVariablesIndex), normalFromBack.density(newValue));
if (flippedVariablesIndex < length - 1) {
    NormalDistribution normal = getNormalDistribution(getValue(currentSample, flippedVariablesIndex + 1));
    double oldForwardPotential = normal.density(oldValue);
    double newForwardPotential = normal.density(newValue);
    // println("Removing old forward potential " + oldForwardPotential);
    currentSample.removePotential(new DoublePotential(oldForwardPotential));
    // println("Multiplying new forward potential " + newForwardPotential);
    currentSample.updatePotential(new DoublePotential(newForwardPotential));
}

// println("Removing old backward potential " + previousP);
currentSample.removePotential(new DoublePotential(previousP));
// println("Multiplying (removing from divisor) old q " + previousQ);
currentSample.updatePotential(new DoublePotential(previousQ));

println("Final sample: " + currentSample);
println();

// check by comparison to brute force calculation of weight:
double productOfPs = 1.0;
for (int i = 1; i != length; i++) {
    productOfPs *= getNormalDistribution(getValue(currentSample, i - 1)).density(getValue(currentSample, i));
}
double productOfQs = Util.fold(fromVariableToQ.values(), (p1, p2) -> p1*p2, 1.0);
double weight = productOfPs/productOfQs;
if (Math.abs(weight - currentSample.getPotential().doubleValue()) > 0.0000001) {
    println("Error in weight calculation");
    System.exit(0);
}

Ini adalah ide yang menarik, tetapi saya melihat beberapa kesulitan dengan itu:

1. bertentangan dengan sampel kepentingan standar, atau bahkan pengambilan sampel kepentingan Metropolised proposal tidak bertindak dalam ruang yang sama dengan distribusi target, tetapi dalam ruang dimensi yang lebih kecil sehingga validasi tidak jelas [dan mungkin memaksakan untuk menjaga bobot lintas iterasi, karenanya menghadapi degenerasi]
2. konstanta normalisasi yang hilang dalam kondisi penuh berubah pada setiap iterasi tetapi tidak diperhitungkan [lihat di bawah]
3. bobot tidak dibatasi, dalam iterasi sepanjang itu, pada akhirnya akan ada simulasi dengan bobot yang sangat besar, kecuali jika seseorang melacak kemunculan terakhir dari pembaruan untuk indeks yang sama , yang mungkin berbenturan dengan validasi Markovian dari sampler Gibbs . Menjalankan percobaan sederhana dengan iterasi dan menunjukkan rentang bobot dari hingga .$j$$n=2$$T=10^3$7.656397e-073.699364e+04

Untuk masuk ke perincian lebih lanjut, pertimbangkan target dua dimensi , termasuk konstanta normalisasi yang tepat, dan terapkan pentingnya Gibbs sampler dengan proposal dan . Bobot kepentingan yang benar [dalam arti menghasilkan ekspektasi yang benar, yaitu, estimator yang tidak bias, untuk fungsi arbitrer dari ] untuk simulasi berturut-turut adalah mana dan adalah margin dari . Atau sederajat $p(\cdot,\cdot)$$q_X(\cdot|y)$$q_Y(\cdot|x)$$(X,Y)$

$$\dfrac{p(x_t,y_{t-1})}{q_X(x_t|y_{t-1})m_Y(y_{t-1})}\qquad\text{or}\qquad\dfrac{p(x_{t-1},y_{t})}{q_Y(y_t|x_{t-1})m_X(x_{t-1})}$$ $m_X(\cdots)$$m_Y(\cdot)$$p(\cdot,\cdot)$ $$\dfrac{p_X(x_t|y_{t-1})}{q_X(x_t|y_{t-1})}\qquad\text{or}\qquad\dfrac{p_Y(y_{t}|x_{t-1})}{q_Y(y_t|x_{t-1})}$$ Dalam kedua kasus ini, ini membutuhkan kepadatan marjinal [tidak bisa diolah] dari dan bawah target .$X$$Y$$p(\cdot,\cdot)$

Penting untuk membandingkan apa yang terjadi di sini dengan algoritma Metropolis tertimbang kepentingan paralel . (Lihat misalnya Schuster und Klebanov, 2018. ) Jika targetnya lagi dan proposal adalah , bobot pentingnya benar [untuk menghasilkan estimasi yang tidak bias] dan tidak memperbarui bobot sebelumnya tetapi mulai dari awal di setiap iterasi.$p(\cdot,\cdot)$$q(\cdot,\cdot|x,y)$

$$\dfrac{p(x',y')}{q(x',y'|x,y)}$$

(C.) Koreksi terhadap kepentingan asli proposal Gibbs adalah mengusulkan nilai baru untuk seluruh vektor, misalnya, , dari proposal Gibbs , karena dengan demikian bobot pentingnya benar [hilang kemungkinan normalisasi konstan yang sekarang benar-benar konstan dan tidak membawa dari iterasi Gibbs sebelumnya] .$(x,y)$$q_X(x_t|y_{t-1})q_Y(y_t|x_{t})$

$$\dfrac{p(x_t,y_t)}{q_X(x_t|y_{t-1})q_Y(y_t|x_{t})}$$

Catatan terakhir: untuk target jalan acak yang dipertimbangkan dalam kode, simulasi langsung dapat dilakukan dengan Cascading: mensimulasikan , kemudian diberikan , & tc.$X_1$$X_2$$X_1$