Menerapkan inferensi variatif stokastik untuk Bayesian Mixture of Gaussian

Saya mencoba untuk menerapkan model Campuran Gaussian dengan inferensi variasional stokastik, berikut ini kertas .

Ini adalah pgm dari Campuran Gaussian.

Menurut makalah itu, algoritma penuh inferensi variatif stokastik adalah:

Dan saya masih sangat bingung dengan metode untuk menskalakannya menjadi GMM.

Pertama, saya pikir parameter variasional lokal hanya dan yang lainnya adalah semua parameter global. Harap perbaiki saya jika saya salah. Apa yang dimaksud dengan langkah 6 ? Apa yang harus saya lakukan untuk mencapai ini? $q_z$ as though Xi is replicated by N times

Bisakah Anda membantu saya dengan ini? Terima kasih sebelumnya!

machine-learning bayesian clustering gaussian-mixture variational-bayes pengguna5779223
sumber

Dikatakan bahwa alih-alih menggunakan seluruh dataset, sampel satu datapoint dan berpura-puralah Anda memiliki datapoint dengan ukuran yang sama. Dalam banyak kasus, ini akan setara dengan mengalikan harapan dengan satu datapoint oleh .

N

$N$

N

$N$

Daeyoung Lim

@ DavideyLim Terima kasih atas balasan Anda! Saya mengerti maksud Anda sekarang, tetapi saya masih bingung bahwa statistik mana yang harus diperbarui secara lokal dan mana yang harus diperbarui secara global. Sebagai contoh, di sini adalah implementasi dari campuran Gaussian, bisakah Anda memberi tahu saya bagaimana skala untuk svi? Saya sedikit tersesat. Terima kasih banyak!

user5779223

Saya tidak membaca seluruh kode tetapi jika Anda berurusan dengan model campuran Gaussian, variabel indikator komponen campuran harus menjadi variabel lokal karena masing-masing terkait dengan hanya satu pengamatan. Jadi variabel laten komponen campuran yang mengikuti distribusi Multinoulli (juga dikenal sebagai distribusi Kategorikal dalam ML) adalah dalam uraian Anda di atas.

z_{i}, i = 1, \dots, N

$z_{i}, \; i=1,\ldots,N$

Daeyoung Lim

@ DavideyLim Ya, saya mengerti apa yang Anda katakan sejauh ini. Jadi untuk distribusi variasional q (Z) q (\ pi, \ mu, \ lambda), q (Z) harus merupakan variabel lokal. Tetapi ada banyak parameter yang terkait dengan q (Z). Di sisi lain, ada juga banyak parameter yang terkait dengan q (\ pi, \ mu, \ lambda). Dan saya tidak tahu cara memperbaruinya dengan tepat.

user5779223

Anda harus menggunakan asumsi bidang-rata untuk mendapatkan distribusi variasi yang optimal untuk parameter variasi. Berikut rujukannya: maths.usyd.edu.au/u/jormerod/JTOpapers/Ormerod10.pdf

Daeyoung Lim

Jawaban:

Tutorial ini ( https://chrisdxie.files.wordpress.com/2016/06/in-depth-variational-inference-tutorial.pdf ) menjawab sebagian besar pertanyaan Anda, dan mungkin akan lebih mudah dipahami daripada makalah SVI asli seperti itu berjalan secara khusus melalui semua detail penerapan SVI (dan mengoordinasikan pendakian VI dan gibbs sampling) untuk model campuran Gaussian (dengan varian yang diketahui).

aleshing
sumber

Pertama, beberapa catatan yang membantu saya memahami makalah SVI:

Dalam menghitung nilai tengah untuk parameter variasional dari parameter global, kami mengambil sampel satu titik data dan berpura-pura seluruh kumpulan data kami ukuran adalah titik tunggal, kali. $N$ $N$
$\eta_g$ adalah parameter alami untuk kondisi penuh variabel global . Notasi digunakan untuk menekankan bahwa ini adalah fungsi dari variabel terkondisi, termasuk data yang diamati. $\beta$

Dalam campuran Gaussians, parameter global kami adalah parameter rata-rata dan presisi (varian terbalik) params untuk masing-masing. Yaitu, adalah parameter alami untuk distribusi ini, sebuah Normal-Gamma dari formulir $k$ $\mu_k, \tau_k$ $\eta_g$

μ, τ \sim N (μ | γ, τ (2 α - 1) G a (τ | α, β)

$\mu, \tau \sim N(\mu|\gamma, \tau(2\alpha -1)Ga(\tau|\alpha, \beta)$

dengan , dan . (Bernardo dan Smith, Bayesian Theory ; perhatikan ini sedikit berbeda dari empat-parameter Normal-Gamma yang biasa Anda lihat .) Kita akan menggunakan untuk merujuk pada parameter variasi untuk $\eta_0 = 2\alpha - 1$ $\eta_1 = \gamma*(2\alpha -1)$ $\eta_2 = 2\beta+\gamma^2(2\alpha-1)$ $a, b, m$ $\alpha, \beta, \mu$

penuh dari adalah Normal-Gamma dengan params , , , di mana adalah yang sebelumnya. (The di sana juga bisa membingungkan; masuk akal dimulai dengan trik diterapkan pada , dan diakhiri dengan jumlah aljabar yang tersisa untuk pembaca.) $\mu_k, \tau_k$ $\dot\eta + \langle\sum_Nz_{n,k}$ $\sum_N z_{n,k}x_N$ $\sum_Nz_{n,k}x^2_{n}\rangle$ $\dot\eta$ $z_{n,k}$ $\exp\ln(p))$ $\prod_N p(x_n|z_n, \alpha, \beta, \gamma) = \prod_N\prod_K\big(p(x_n|\alpha_k,\beta_k,\gamma_k)\big)^{z_{n,k}}$

Dengan itu, kita dapat menyelesaikan langkah (5) dari pseudocode SVI dengan:

ϕ_{n, k} \propto \exp (l n (π) + E_{q} \ln (p (x_{n} | α_{k}, β_{k}, γ_{k})) = \exp (\ln (π) + E_{q} [⟨ μ_{k} τ_{k}, \frac{- τ}{2} ⟩ \cdot ⟨ x, x^{2} ⟩ - \frac{μ^{2} τ - \ln τ}{2})]

$\phi_{n,k} \propto \exp (ln(\pi) + \mathbb E_q \ln(p(x_n|\alpha_k, \beta_k, \gamma_k))\\ =\exp(\ln(\pi) + \mathbb E_q \big[\langle \mu_k\tau_k, \frac{-\tau}{2} \rangle \cdot\langle x, x^2\rangle - \frac{\mu^2\tau - \ln \tau}{2})\big]$

Memperbarui parameter global lebih mudah, karena setiap parameter terkait dengan jumlah data atau salah satu statistik yang memadai:

\hat{λ} = \dot{η} + N ϕ_{n} ⟨ 1, x, x^{2} ⟩

$\hat \lambda = \dot \eta + N\phi_n \langle 1, x, x^2 \rangle$

Inilah kemungkinan marginal dari data terlihat pada banyak iterasi, ketika dilatih tentang data yang sangat tiruan, mudah dipisahkan (kode di bawah). Plot pertama menunjukkan kemungkinan dengan inisialisasi, parameter variasional acak dan iterasi; masing-masing berikutnya adalah setelah kekuatan dua iterasi berikutnya. Dalam kode, merujuk ke parameter variasi untuk . $0$ $a, b, m$ $\alpha, \beta, \mu$

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Sun Aug 12 12:49:15 2018

@author: SeanEaster
"""

import numpy as np
from matplotlib import pylab as plt
from scipy.stats import t
from scipy.special import digamma 

# These are priors for mu, alpha and beta

def calc_rho(t, delay=16,forgetting=1.):
    return np.power(t + delay, -forgetting)

m_prior, alpha_prior, beta_prior = 0., 1., 1.
eta_0 = 2 * alpha_prior - 1
eta_1 = m_prior * (2 * alpha_prior - 1)
eta_2 = 2 *  beta_prior + np.power(m_prior, 2.) * (2 * alpha_prior - 1)

k = 3

eta_shape = (k,3)
eta_prior = np.ones(eta_shape)
eta_prior[:,0] = eta_0
eta_prior[:,1] = eta_1
eta_prior[:,2] = eta_2

np.random.seed(123) 
size = 1000
dummy_data = np.concatenate((
        np.random.normal(-1., scale=.25, size=size),
        np.random.normal(0.,  scale=.25,size=size),
        np.random.normal(1., scale=.25, size=size)
        ))
N = len(dummy_data)
S = 1

# randomly init global params
alpha = np.random.gamma(3., scale=1./3., size=k)
m = np.random.normal(scale=1, size=k)
beta = np.random.gamma(3., scale=1./3., size=k)

eta = np.zeros(eta_shape)
eta[:,0] = 2 * alpha - 1
eta[:,1] = m * eta[:,0]
eta[:,2] = 2. * beta + np.power(m, 2.) * eta[:,0]


phi = np.random.dirichlet(np.ones(k) / k, size = dummy_data.shape[0])

nrows, ncols = 4, 5
total_plots = nrows * ncols
total_iters = np.power(2, total_plots - 1)
iter_idx = 0

x = np.linspace(dummy_data.min(), dummy_data.max(), num=200)

while iter_idx < total_iters:

    if np.log2(iter_idx + 1) % 1 == 0:

        alpha = 0.5 * (eta[:,0] + 1)
        beta = 0.5 * (eta[:,2] - np.power(eta[:,1], 2.) / eta[:,0])
        m = eta[:,1] / eta[:,0]
        idx = int(np.log2(iter_idx + 1)) + 1

        f = plt.subplot(nrows, ncols, idx)
        s = np.zeros(x.shape)
        for _ in range(k):
            y = t.pdf(x, alpha[_], m[_], 2 * beta[_] / (2 * alpha[_] - 1))
            s += y
            plt.plot(x, y)
        plt.plot(x, s)
        f.axes.get_xaxis().set_visible(False)
        f.axes.get_yaxis().set_visible(False)

    # randomly sample data point, update parameters
    interm_eta = np.zeros(eta_shape)
    for _ in range(S):
        datum = np.random.choice(dummy_data, 1)

        # mean params for ease of calculating expectations
        alpha = 0.5 * ( eta[:,0] + 1)
        beta = 0.5 * (eta[:,2] - np.power(eta[:,1], 2) / eta[:,0])
        m = eta[:,1] / eta[:,0]

        exp_mu = m
        exp_tau = alpha / beta 
        exp_tau_m_sq = 1. / (2 * alpha - 1) + np.power(m, 2.) * alpha / beta
        exp_log_tau = digamma(alpha) - np.log(beta)


        like_term = datum * (exp_mu * exp_tau) - np.power(datum, 2.) * exp_tau / 2 \
            - (0.5 * exp_tau_m_sq - 0.5 * exp_log_tau)
        log_phi = np.log(1. / k) + like_term
        phi = np.exp(log_phi)
        phi = phi / phi.sum()

        interm_eta[:, 0] += phi
        interm_eta[:, 1] += phi * datum
        interm_eta[:, 2] += phi * np.power(datum, 2.)

    interm_eta = interm_eta * N / S
    interm_eta += eta_prior

    rho = calc_rho(iter_idx + 1)

    eta = (1 - rho) * eta + rho * interm_eta

    iter_idx += 1

Sean Easter
sumber