spm_vb_get_Gn

 Compute Gn for VB GLMAR modelling 
 FORMAT [G,G1,G2,G3] = spm_vb_get_Gn (Y,slice,n)

 Y             [T x N] time series 
 slice         data structure 
 n             voxel number

 %W% Will Penny and Nelson Trujillo-Barreto %E%

0001 function [G,G1,G2,G3] = spm_vb_get_Gn (Y,slice,n)
0002 % Compute Gn for VB GLMAR modelling
0003 % FORMAT [G,G1,G2,G3] = spm_vb_get_Gn (Y,slice,n)
0004 %
0005 % Y             [T x N] time series
0006 % slice         data structure
0007 % n             voxel number
0008 %
0009 % %W% Will Penny and Nelson Trujillo-Barreto %E%
0010 
0011 p=slice.p;
0012 k=slice.k;
0013 
0014 a_mean_t=slice.ap_mean(:,n);
0015 a_mean=a_mean_t';
0016 block_n = [(n-1)*k+1:n*k];
0017 w_mean=slice.w_mean(block_n,1);
0018 w2=slice.w2{n};
0019 a2=slice.a2{n};
0020 
0021 % Equations 77 to 81 of paper VB1 but implemented
0022 % efficiently using cross-covariance method described in paper VB3
0023 G11=slice.y2(n);
0024 G12=ones(1,p)*(a2.*slice.I.Gy(:,:,n))*ones(p,1);
0025 G13=-2*slice.I.gy(:,n)'*a_mean_t;
0026 G1=G11+G12+G13;
0027 
0028 G21=ones(1,k)*(w2.*slice.I.Gx)*ones(k,1);
0029 G22=trace(slice.I.A2_tilde(:,:,n)*slice.w_cov{n});
0030 G23=ones(1,k)*((w_mean*w_mean').*slice.I.A2_tilde(:,:,n))*ones(k,1);
0031 G24=2*ones(1,k)*(w2.*slice.I.A3a_tilde(:,:,n))*ones(k,1);
0032 G2=G21+G22+G23+G24;
0033 
0034 G31=-2*w_mean'*slice.I.gxy(:,n);
0035 G32=2*a_mean*slice.I.rxy(:,:,n)*w_mean;
0036 G33=2*a_mean*slice.I.Gxy(:,:,n)*w_mean;
0037 G34=-2*ones(1,p)*(a2.*slice.I.W_tilde(:,:,n))*ones(p,1);
0038 G3=G31+G32+G33+G34;
0039 
0040 G=G1+G2+G3;
0041 
0042