fcalc.m

基于变分贝叶斯方法的状态空间模型学习程序

% Function to calculate the lower bound F, for the variational Bayesian state-space model.
%
% Broken down into 6 parts.  These should be obvious with reference to thesis.
%
% Matthew J. Beal GCNU 01/03/03

function [F] = Fcalc(pa,pb,alpha,beta,gamma,delta,SigA,SigB,SigC,SigD,Bb,Db,SA,SC,GA,MD,G,lnZpn,Tn)

F = zeros(1,6);
k = size(alpha,1);
pinp = size(beta,1);
p = size(G,1);

T = sum(Tn,2); 

F(1) = -( -k/2*( log(prod(beta))+lndet(SigB) ) ...
	  + k/2*trace( diag(beta)*SigB-eye(pinp,pinp) ) ...
	  + 1/2*trace( diag(beta)*SigB*Bb*Bb'*SigB ) );
F(2) = -( -k/2*( log(prod(alpha))+lndet(SigA) )...
	  + k/2*trace( diag(alpha)*SigA-eye(k,k) ) ...
	  + k/2*trace( diag(alpha)*SigA*GA*SigB*GA'*SigA ) ...
	  +1/2*trace( diag(alpha)*SigA*(SA-GA*SigB*Bb)*(SA-GA*SigB*Bb)'*SigA ) );
F(3) = - klgamma(pa+sum(Tn,2)/2,pb+G'/2,pa,pb);
F(4) = -( -p/2*( log(prod(delta))+lndet(SigD) ) ...
	  + p/2*trace( diag(delta)*SigD-eye(pinp,pinp) ) ...
	  + 1/2*trace( diag(delta)*SigD*Db*(pa+sum(Tn,2)/2)*diag(1./(pb+G/2))*Db'*SigD ) );
F(5) = -( -p/2*( log(prod(gamma))+lndet(SigC) ) ...
	  + p/2*trace( diag(gamma)*SigC-eye(k,k) ) ...
	  + p/2*trace( diag(gamma)*SigC*MD*SigD*MD'*SigC ) ...
	  + 1/2*trace( diag(gamma)*SigC*(SC-MD*SigD*Db)*(pa+sum(Tn,2)/2)*diag(1./(pb+G/2))*(SC-MD*SigD*Db)'*SigC ) );
F(6) = sum(lnZpn,2);