/
spm_kl_normald.m
40 lines (35 loc) · 1.19 KB
/
spm_kl_normald.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
function [d] = spm_kl_normald (m_q,c_q,m_p,c_p)
% KL divergence between two Gaussians with, possibly, diagonal covariances
% FORMAT [d] = spm_kl_normald (m_q,c_q,m_p,c_p)
%
% Calculate the KL distance KL (Q||P) = <log Q/P> where avg is wrt Q
% between two Normal densities Q and P
%
% m_q, c_q Mean and covariance of first Normal density
% m_p, c_p Mean and covariance of second Normal density
%
% If c_q and c_p are diagonal, pass them as vectors, and KL will
% be computed more efficiently. Both must be full or both must be diagonal
%___________________________________________________________________________
% Copyright (C) 2008 Wellcome Trust Centre for Neuroimaging
% Will Penny
% $Id: spm_kl_normald.m 2696 2009-02-05 20:29:48Z guillaume $
dd=length(m_q);
m_q=m_q(:);
m_p=m_p(:);
[d1,d2]=size(c_q);
if (d1==d2)
% Full covariances
Term1=0.5*spm_logdet(c_p)-0.5*spm_logdet(c_q);
inv_c_p=inv(c_p);
Term2=0.5*trace(inv_c_p*c_q)+0.5*(m_q-m_p)'*inv_c_p*(m_q-m_p);
else
% Diagonal covariances
c_q=c_q(:);
c_p=c_p(:);
Term1=0.5*sum(log(c_p))-0.5*sum(log(c_q));
Term2=0.5*sum(c_q./c_p);
e=sqrt(1./c_p).*(m_q-m_p);
Term2=Term2+0.5*e'*e;
end
d=Term1+Term2-0.5*dd;