变分推断(Variational Inference)解析

一、什么是变分推断

假设在一个贝叶斯模型中, x x x为一组观测变量, z z z为一组隐变量(参数也看做随机变量,包含在 z z z中),则推断问题为计算后验概率密度 P = ( z ∣ x ) P=(z|x) P=(zx)。根据贝叶斯公式,有:
p ( z ∣ x ) = p ( x , z ) p ( x ) = p ( x , z ) ∫ p ( x , z ) d z p(z|x)=\frac{p(x,z)}{p(x)}=\frac{p(x,z)}{\int p(x,z)dz} p(zx)=p(x)p(x,z)=p(x,z)dzp(x,z)
但是在实际应用中,可能由于积分没有闭式解,或者是指数级的计算复杂度等原因,导致计算上面公式中的积分往往是不可行的。变分推断就是用来解决这个问题的。

变分推断是变分法在推断问题中的应用,既然无法直接求得后验概率密度 p ( z ∣ x ) p(z|x) p(zx),那我们可以寻找一个简单的分布 q ∗ ( z ) q^*(z) q(z)来近似后验概率密度 p ( z ∣ x ) p(z|x) p(zx),这就是变分推断的思想。借此,我们将推断问题转换为一个泛函优化问题:
q ∗ ( z ) = arg ⁡ min ⁡ q ( z ) ∈ Q K L ( q ( z ) ∣ ∣ p ( z ∣ x ) ) (1) q^*(z)=\arg\min_{q(z)\in Q}KL(q(z)||p(z|x))\tag{1} q(z)=argq(z)QminKL(q(z)p(zx))(1)
其中 Q Q Q为候选的概率分布族。但是又出现了一个新的问题:我们已经知道后验概率密度 p ( z ∣ x ) p(z|x) p(zx)难以计算,所以上式中的KL散度本身也是无法计算的!这时,需要借助于证据下界ELBO。

ELBO

ELBO,全称为 Evidence Lower Bound,即证据下界。这里的证据指数据或可观测变量的概率密度。

假设 x = x 1 : n x=x_{1:n} x=x1:n表示一系列可观测数据集, z = z 1 : m z=z_{1:m} z=z1:m为一系列隐变量(latent variables)。则可用 p ( z , x ) p(z,x) p(z,x)表示联合概率, p ( z ∣ x ) p(z∣x) p(zx)为条件概率, p ( x ) p(x) p(x)为证据。

那么,贝叶斯推理需要求解的就是条件概率,即: p ( z ∣ x ) = p ( x , z ) p ( x ) p(z|x)=\frac{p(x,z)}{p(x)} p(zx)=p(x)p(x,z)
(1)式中的KL散度可以表示为 K L ( q ( z ) ∣ ∣ p ( z ∣ x ) ) = ∫ q ( z ) log ⁡ q ( z ) p ( z ∣ x ) d z KL(q(z)||p(z|x))=\int q(z)\log\frac{q(z)}{p(z|x)}dz KL(q(z)p(zx))=q(z)logp(zx)q(z)dz其中, x x x为可观测数据集, z z z为未知变量,下面将公式继续变形:
∫ q ( z ) log ⁡ q ( z ) p ( z ∣ x ) d z = − ∫ q ( z ) log ⁡ p ( z ∣ x ) q ( z ) d z = − ∫ q ( z ) log ⁡ p ( x , z ) q ( z ) p ( x ) d z = − ∫ q ( z ) log ⁡ p ( x , z ) d z + ∫ q ( z ) log ⁡ q ( z ) d z + ∫ q ( z ) log ⁡ p ( x ) d z \begin{aligned}\int q(z)\log\frac{q(z)}{p(z|x)}dz&=-\int q(z)\log\frac{p(z|x)}{q(z)}dz\\&=-\int q(z)\log\frac{p(x,z)}{q(z)p(x)}dz\\&=-\int q(z)\log p(x,z)dz+\int q(z)\log q(z)dz+\int q(z)\log p(x)dz\end{aligned} q(z)logp(zx)q(z)dz=q(z)logq(z)p(zx)dz=q(z)logq(z)p(x)p(x,z)dz=q(z)logp(x,z)dz+q(z)logq(z)dz+q(z)logp(x)dz其中, ∫ q ( z ) d z = 1 \int q(z)dz=1 q(z)dz=1进而可以转化成: = − ∫ q ( z ) log ⁡ p ( x , z ) d z + ∫ q ( z ) log ⁡ q ( z ) d z + log ⁡ p ( x ) =-\int q(z)\log p(x,z)dz+\int q(z)\log q(z)dz+\log p(x) =q(z)logp(x,z)dz+q(z)logq(z)dz+logp(x) L ( q ( z ) ) = ∫ q ( z ) log ⁡ p ( x , z ) d z − ∫ q ( z ) log ⁡ q ( z ) d z L(q(z))=\int q(z)\log p(x,z)dz-\int q(z)\log q(z)dz L(q(z))=q(z)logp(x,z)dzq(z)logq(z)dz
则有 K L ( q ( z ) ∣ ∣ p ( z ∣ x ) ) = − L ( q ( z ) ) + log ⁡ p ( x ) KL(q(z)||p(z|x))=-L(q(z))+\log p(x) KL(q(z)p(zx))=L(q(z))+logp(x)从这个公式可以发现, log ⁡ p ( x ) \log p(x) logp(x)不涉及参数(数据似然),因此在最小化 K L ( q ( z ) ∣ ∣ p ( z ∣ x ) ) KL(q(z)||p(z|x)) KL(q(z)p(zx))时可以忽略。那么,最小化 K L ( q ( z ) ∣ ∣ p ( z ∣ x ) ) KL(q(z)||p(z|x)) KL(q(z)p(zx))便转化成了最大化 L ( q ( z ) ) L(q(z)) L(q(z))

因为 K L ( q ( z ) ∣ ∣ p ( z ∣ x ) ) ≥ 0 KL(q(z)||p(z|x))\geq 0 KL(q(z)p(zx))0,即: − L ( q ( z ) ) + log ⁡ p ( x ) ≥ 0 -L(q(z))+\log p(x)\geq 0 L(q(z))+logp(x)0进而可以得到: log ⁡ p ( x ) ≥ L ( q ( z ) ) \log p(x)\geq L(q(z)) logp(x)L(q(z))因此,可以将 L ( q ( z ) ) L(q(z)) L(q(z))堪称 log ⁡ p ( x ) \log p(x) logp(x)的下界,这个下界也称之为ELBO(evidence lower bound),那么最小化 K L ( q ( z ) ∣ ∣ p ( z ∣ x ) ) KL(q(z)||p(z|x)) KL(q(z)p(zx)),可以看成最大化下界的问题。

另外,从公式中可以看到,KL散度是 L ( q ( z ) ) L(q(z)) L(q(z)) log ⁡ p ( x ) \log p(x) logp(x)的误差,当然误差越小越好。

根据以上结果,最新的目标函数转化成了 q ∗ ( z ) = arg ⁡ max ⁡ q ( z ) ∈ Q L ( q ( z ) ) = arg ⁡ max ⁡ q ( z ) ∈ Q ∫ z q ( z ) log ⁡ p ( x , z ) d z ⏟ ( a ) − ∫ z q ( z ) log ⁡ q ( z ) d z ⏟ ( b ) (2) \begin{aligned}q^*(z)&=\arg\max_{q(z)\in Q}L(q(z))\\&=\arg\max_{q(z)\in Q}\underbrace{\int_z q(z)\log p(x,z)dz}_{(a)}-\underbrace{\int_z q(z)\log q(z)dz}_{(b)}\tag{2}\end{aligned} q(z)=argq(z)QmaxL(q(z))=argq(z)Qmax(a) zq(z)logp(x,z)dz(b) zq(z)logq(z)dz(2)至此,我们已经解决了KL散度无法求解的问题,将泛函优化问题转换为寻找一个简单分布 q ∗ ( z ) q^*(z) q(z)来最大化证据下界 L ( q ( z ) ) L(q(z)) L(q(z))

二、基于平均场理论的变分推断

在变分推断中,候选分布族 Q Q Q的复杂性决定了优化问题的复杂性。一个通常的选择是平均场分布族,即 z z z可以拆分成多组相互独立的变量,有: q ( z ) = ∏ i = 1 M q i ( z i ) (3) q(z)=\prod_{i=1}^Mq_i(z_i)\tag{3} q(z)=i=1Mqi(zi)(3)其中 z i z_i zi是隐变量的子集,可以是单变量,也可以是一组多元变量。

下面我们分布(3)把将代入(2)中的(a)和(b),看看 L ( q ( z ) ) L(q(z)) L(q(z))最后的模样,其中假设我们想先求 q j ( z j ) q_j(z_j) qj(zj),将其它组的 q ∖ j ( z ∖ j ) q_{\setminus j}(z_{\setminus j}) qj(zj)当作常量:

2.1、求解(a)

我们首先求解(a):
( a ) = ∫ z q ( z )

  • 15
    点赞
  • 58
    收藏
    觉得还不错? 一键收藏
  • 3
    评论
作为一个AI语言模型,我并不具备编写Matlab程序的能力。但是我可以为您提供一些变分推断的Matlab程序示例,供您参考。 1. 高斯混合模型的变分推断 %生成数据 n = 1000; % 数据量 k = 3; % 高斯分量个数 d = 2; % 数据维度 % 生成高斯分量的均值和协方差矩阵 mu = rand(k,d)*10-5; % 均值矩阵 sigma = zeros(d,d,k); % 协方差矩阵 for i=1:k sigma(:,:,i) = diag(rand(d,1)*5); % 对角线矩阵 end % 生成数据 data = zeros(n,d); for i=1:n k_index = randi(k); data(i,:) = mvnrnd(mu(k_index,:),sigma(:,:,k_index)); end % 变分推断 % 初始化参数 alpha = ones(1,k)./k; % 高斯分量的先验概率 mu_init = rand(k,d)*10-5; % 均值矩阵 sigma_init = zeros(d,d,k); % 协方差矩阵 for i=1:k sigma_init(:,:,i) = diag(rand(d,1)*5); % 对角线矩阵 end q_alpha = alpha; % 变分分布的先验概率 q_mu = mu_init; % 均值矩阵 q_sigma = sigma_init; % 协方差矩阵 % 迭代计算 max_iter = 100; for iter=1:max_iter % 更新q_alpha E_z = zeros(n,k); for i=1:n for j=1:k E_z(i,j) = log(alpha(j)) + log_mvnpdf(data(i,:),q_mu(j,:),q_sigma(:,:,j)); end E_z(i,:) = exp(E_z(i,:) - max(E_z(i,:))); % 防止指数爆炸 E_z(i,:) = E_z(i,:) ./ sum(E_z(i,:)); end q_alpha = alpha + sum(E_z,1); % 更新q_mu for j=1:k q_mu(j,:) = sum(E_z(:,j).*data,1) ./ sum(E_z(:,j)); end % 更新q_sigma for j=1:k diff = data - q_mu(j,:); q_sigma(:,:,j) = diff'*diag(E_z(:,j))*diff ./ sum(E_z(:,j)); end end % 计算后验概率 posterior = zeros(n,k); for i=1:n for j=1:k posterior(i,j) = log(q_alpha(j)) + log_mvnpdf(data(i,:),q_mu(j,:),q_sigma(:,:,j)); end posterior(i,:) = exp(posterior(i,:) - max(posterior(i,:))); % 防止指数爆炸 posterior(i,:) = posterior(i,:) ./ sum(posterior(i,:)); end % 显示结果 figure; hold on; scatter(data(:,1),data(:,2),10,posterior(:,1),'filled'); scatter(q_mu(:,1),q_mu(:,2),100,'k','filled'); scatter(mu(:,1),mu(:,2),100,'r','filled'); hold off; title('GMM with Variational Inference'); legend('Cluster 1','Cluster 2','Cluster 3'); xlabel('Feature 1'); ylabel('Feature 2'); 2. 隐马尔可夫模型的变分推断 % 生成数据 n = 1000; % 数据量 k = 3; % 隐状态个数 d = 2; % 数据维度 % 生成隐状态转移矩阵和观测矩阵 A = rand(k,k); % 隐状态转移矩阵 A = A ./ sum(A,2); B = rand(k,d)*10-5; % 观测矩阵 % 生成数据 data = zeros(n,d); z = zeros(n,1); z(1) = randi(k); data(1,:) = mvnrnd(B(z(1),:),eye(d)); for i=2:n z(i) = randsample(k,1,true,A(z(i-1),:)); data(i,:) = mvnrnd(B(z(i),:),eye(d)); end % 变分推断 % 初始化参数 alpha = ones(1,k)./k; % 隐状态的先验概率 A_init = rand(k,k); % 隐状态转移矩阵 A_init = A_init ./ sum(A_init,2); B_init = rand(k,d)*10-5; % 观测矩阵 q_alpha = alpha; % 变分分布的先验概率 q_A = A_init; % 隐状态转移矩阵 q_B = B_init; % 观测矩阵 % 迭代计算 max_iter = 100; for iter=1:max_iter % 更新q_alpha E_z = zeros(n,k); for i=1:n for j=1:k E_z(i,j) = log(alpha(j)) + log(A(z(i-1),j)) + log_mvnpdf(data(i,:),q_B(j,:),eye(d)); end E_z(i,:) = exp(E_z(i,:) - max(E_z(i,:))); % 防止指数爆炸 E_z(i,:) = E_z(i,:) ./ sum(E_z(i,:)); end q_alpha = alpha + E_z(1,:); % 更新q_A for i=1:k for j=1:k q_A(i,j) = sum(E_z(1:end-1,i).*E_z(2:end,j)) ./ sum(E_z(1:end-1,i)); end end % 更新q_B for j=1:k q_B(j,:) = sum(E_z(:,j).*data,1) ./ sum(E_z(:,j)); end end % 计算后验概率 posterior = zeros(n,k); for i=1:n for j=1:k if i==1 posterior(i,j) = log(q_alpha(j)) + log_mvnpdf(data(i,:),q_B(j,:),eye(d)); else posterior(i,j) = log(A(z(i-1),j)) + log_mvnpdf(data(i,:),q_B(j,:),eye(d)); end end posterior(i,:) = exp(posterior(i,:) - max(posterior(i,:))); % 防止指数爆炸 posterior(i,:) = posterior(i,:) ./ sum(posterior(i,:)); end % 显示结果 figure; hold on; scatter(data(:,1),data(:,2),10,posterior(:,1),'filled'); scatter(q_B(:,1),q_B(:,2),100,'k','filled'); scatter(B(:,1),B(:,2),100,'r','filled'); hold off; title('HMM with Variational Inference'); legend('State 1','State 2','State 3'); xlabel('Feature 1'); ylabel('Feature 2');

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 3
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值