| 𝐙={Z1…Zn} given some data P(𝐙∣𝐗)≈Q(𝐙).
The distribution Q(𝐙) similar to the true posterior, Q(𝐙) that minimizes P(𝐙∣𝐗)=P(𝐗∣𝐙)P(𝐙)P(𝐗)=P(𝐗∣𝐙)P(𝐙)∫𝐙P(𝐗,𝐙′)d𝐙′
The marginalization over 𝐙 is combinatorially large. Therefore, we seek an approximation, using
Given that P(𝐗) is a constant with respect to Q(𝐙) is a distribution, we have
- logP(𝐗)=DKL(Q∥P)−𝔼𝐐[logQ(𝐙)−logP(𝐙,𝐗)]=DKL(Q∥P)+ℒ(Q)
As the log-evidence ℒ(Q) minimizes the KL divergence of Q, P(𝐙∣𝐗), and a lower bound ℒ(Q) is known as the (negative) variational free energy in analogy with thermodynamic free energy because it can also be expressed as a negative energy ℒ(Q) is also known as Evidence Lower Bound, abbreviated as ELBO, to emphasize that it is a lower (worst-case) bound on the log-evidence of the data.
By the generalized Pythagorean theorem of Bregman divergence, of which KL-divergence is a special case, it can be shown that:
 Generalized Pythagorean theorem for Bregman divergence
- Q=Q∗≜argminQ∈𝒞DKL(Q∥P).
In this case, the global minimizer q∗(𝐙2)=P(𝐗)ζ(𝐗)P(𝐙2∣𝐗)exp(DKL(q∗(𝐙1∣𝐙2)∥P(𝐙1∣𝐙2,𝐗)))=1ζ(𝐗)exp𝔼q∗(𝐙1∣𝐙2)(logP(𝐙,𝐗)q∗(𝐙1∣𝐙2)),
in which the normalizing constant is:
- P(𝐗)≥ζ(𝐗)=exp(ℒ(Q∗)), as shown above.
By interchanging the roles of q∗(𝐙1) and P(𝐙2∣𝐗), respectively. Although this iterative scheme is guaranteed to converge monotonically, the converged 𝒞 is confined within independent space, i.e. Q(𝐙) is usually assumed to factorize over some partition of the latent variables, i.e. for some partition of the latent variables Q(𝐙)=∏i=1Mqi(𝐙i∣𝐗)
It can be shown using the calculus of variations (hence the name "variational Bayes") that the "best" distribution qj∗(𝐙j∣𝐗)=eEq−j∗[lnp(𝐙,𝐗)]∫eEq−j∗[lnp(𝐙,𝐗)]d𝐙j
where qj∗(𝐙j∣𝐗).
In practice, we usually work in terms of logarithms, i.e.:
- Eq−j∗[lnp(𝐙,𝐗)] can usually be simplified into a function of the fixed hyperparameters of the prior distributions over the latent variables and of expectations (and sometimes higher moments such as the variance) of latent variables not in the current partition (i.e. latent variables not included in (Θ,ℱ,P) and λ such that h denote any real-valued random variable on logEP[exph]=supQ≪P{EQ[h]−DKL(Q∥P)}.
Further, the supremum on the right-hand side is attained if and only if it holds
- p(θ)=dP/dλ and Q with respect to τ, i.e. the mean also follows a Gaussian distribution while the precision follows a gamma distribution. In other words:
- b0 in the prior distributions are fixed, given values. They can be set to small positive numbers to give broad prior distributions indicating ignorance about the prior distributions of N data points μ and p(𝐗∣μ,τ)=∏n=1N𝒩(xn∣μ,τ−1)p(μ∣τ)=𝒩(μ∣μ0,(λ0τ)−1)p(τ)=Gamma(τ∣a0,b0)
where
- μ and C, μ. Note that the term μ. Hence in line 3 we can absorb it into the constant term at the end. We do the same thing in line 7.
The last line is simply a quadratic polynomial in qμ∗(μ) itself is a Gaussian distribution.
With a certain amount of tedious math (expanding the squares inside of the braces, separating out and grouping the terms involving μ), we can derive the parameters of the Gaussian distribution:
- q(τ)edit
The derivation of qτ∗(τ) is a gamma distribution. Specifically:
- qτ∗(τ)∼Gamma(τ∣aN,bN)aN=a0+N+12bN=b0+12Eμ[∑n=1N(xn−μ)2+λ0(μ−μ0)2]
In each case, the parameters for the distribution over one of the variables depend on expectations taken with respect to the other variable. We can expand the expectations, using the standard formulas for the expectations of moments of the Gaussian and gamma distributions:
- bN=b0+12Eμ[∑n=1N(xn−μ)2+λ0(μ−μ0)2]=b0+12Eμ[(λ0+N)μ2−2(λ0μ0+∑n=1Nxn)μ+(∑n=1Nxn2)+λ0μ02]=b0+12[(λ0+N)Eμ[μ2]−2(λ0μ0+∑n=1Nxn)Eμ[μ]+(∑n=1Nxn2)+λ0μ02]=b0+12[(λ0+N)(λN−1+μN2)−2(λ0μ0+∑n=1Nxn)μN+(∑n=1Nxn2)+λ0μ02]
We can then write the parameter equations as follows, without any expectations:
- bN. This naturally suggests an EM-like algorithm:
- Compute μN and λN, along with the known values of the other parameters, to compute λN.
- Repeat the last two steps until convergence (i.e. until neither value has changed more than some small amount).
We then have values for the hyperparameters of the approximating distributions of the posterior parameters, which we can use to compute any properties we want of the posterior — e.g. its mean and variance, a 95% highest-density region (the smallest interval that includes 95% of the total probability), etc.
It can be shown that this algorithm is guaranteed to converge to a local maximum.
Note also that the posterior distributions have the same form as the corresponding prior distributions. We did not assume this; the only assumption we made was that the distributions factorize, and the form of the distributions followed naturally. It turns out (see below) that the fact that the posterior distributions have the same form as the prior distributions is not a coincidence, but a general result whenever the prior distributions are members of the exponential family, which is the case for most of the standard distributions.
The above example shows the method by which the variational-Bayesian approximation to a posterior probability density in a given Bayesian network is derived:
- Describe the network with a graphical model, identifying the observed variables (data) 𝐙) and their conditional probability distributions. Variational Bayes will then construct an approximation to the posterior probability 𝐙j, write down the formula for the best approximating distribution
A more complex exampleedit Bayesian Gaussian mixture model using plate notation. Smaller squares indicate fixed parameters; larger circles indicate random variables. Filled-in shapes indicate known values. The indication [K] means a vector of size K; [D,D] means a matrix of size D×D; K alone means a categorical variable with K outcomes. The squiggly line coming from z ending in a crossbar indicates a switch — the value of this variable selects, for the other incoming variables, which value to use out of the size-K array of possible values.
Imagine a Bayesian Gaussian mixture model described as follows:
- α0. The Dirichlet distribution is the conjugate prior of the categorical distribution or multinomial distribution.
- 𝒩() is the Gaussian distribution, in this case specifically the multivariate Gaussian distribution.
The interpretation of the above variables is as follows:
- D-dimensional vector distributed according to a multivariate Gaussian distribution.
- k=1…K, as described above.
- 𝝁i=1…K and p(𝐗∣𝐙,𝝁,𝜦)=∏n=1N∏k=1K𝒩(𝐱n∣𝝁k,𝜦k−1)znkp(𝐙∣𝝅)=∏n=1N∏k=1Kπkznkp(𝝅)=Γ(Kα0)Γ(α0)K∏k=1Kπkα0−1p(𝝁∣𝜦)=∏k=1K𝒩(𝝁k∣𝝁0,(β0𝜦k)−1)p(𝜦)=∏k=1K𝒲(𝜦k∣𝐖0,ν0)
where
- lnq∗(𝐙)=E𝝅,𝝁,𝜦[lnp(𝐗,𝐙,𝝅,𝝁,𝜦)]+constant=E𝝅[lnp(𝐙∣𝝅)]+E𝝁,𝜦[lnp(𝐗∣𝐙,𝝁,𝜦)]+constant=∑n=1N∑k=1Kznklnρnk+constant
where we have defined
- q∗(𝐙)∝∏n=1N∏k=1Kρnkznk
Requiring that this be normalized ends up requiring that the q∗(𝐙)=∏n=1N∏k=1Krnkznk
where
- 𝐳n, which is distributed as a single-observation multinomial distribution with parameters E[znk]=rnk
which is a standard result for categorical distributions.
Now, considering the factor lnq∗(𝝅)=lnp(𝝅)+E𝐙[lnp(𝐙∣𝝅)]+constant=(α0−1)∑k=1Klnπk+∑n=1N∑k=1Krnklnπk+constant
Taking the exponential of both sides, we recognize αk=α0+Nk
where
- 𝝁k and βk=β0+Nk𝐦k=1βk(β0𝝁0+Nk𝐱¯k)𝐖k−1=𝐖0−1+Nk𝐒k+β0Nkβ0+Nk(𝐱¯k−𝝁0)(𝐱¯k−𝝁0)Tνk=ν0+NkNk=∑n=1Nrnk𝐱¯k=1Nk∑n=1Nrnk𝐱n𝐒k=1Nk∑n=1Nrnk(𝐱n−𝐱¯k)(𝐱n−𝐱¯k)T
Finally, notice that these functions require the values of E[lnπk], E𝝁k,𝜦k[(𝐱n−𝝁k)T𝜦k(𝐱n−𝝁k)]=Dβk−1+νk(𝐱n−𝐦k)T𝐖k(𝐱n−𝐦k)lnΛ~k≡E[ln|𝜦k|]=∑i=1Dψ(νk+1−i2)+Dln2+ln|𝐖k|lnπ~k≡E[ln|πk|]=ψ(αk)−ψ(∑i=1Kαi)
These results lead to
- βk, νk of the variables Nk, rnk.
- The update equations for the parameters Nk, which depends in turn on βk, νk as well as an indirect circular dependence on α1…K through rnk using the current values of all the other parameters.
- An M-step that uses the new value of p(𝐙∣𝐗); the computation of the statistics 𝐒k corresponds closely to the computation of corresponding "soft-count" statistics over the data; and the use of those statistics to compute new values of the parameters corresponds closely to the use of soft counts to compute new parameter values in normal EM over a Gaussian mixture model.
Exponential-family distributionseditNote that in the previous example, once the distribution over unobserved variables was assumed to factorize into distributions over the "parameters" and distributions over the "latent data", the derived "best" distribution for each variable was in the same family as the corresponding prior distribution over the variable. This is a general result that holds true for all prior distributions derived from the exponential family.
|
Generalized Pythagorean theorem for Bregman divergence
