---
function [Vindx, J ,emin] = vb_dipole_estimate(bexp, basis, L, Cov)
単一ダイポール推定
Vindx : 推定されたダイポールのインデックス
: ダイポールの位置座標 = V(Vindx,:)
J : ダイポール電流強度
emin : 復元誤差
bexp : MEG-data (Nsensor ,Time)
basis : Lead Field (Ndipole, Nsensor)
L = 2 (電流方向)
Cov : Noise covariance matrix : Cov(sensor,sensor)
---
Copyright (C) 2011, ATR All Rights Reserved.
License : New BSD License(see VBMEG_LICENSE.txt)0001 function [Vindx, J ,emin] = vb_dipole_estimate_cov(bexp, basis, L, Cov) 0002 % --- 0003 % function [Vindx, J ,emin] = vb_dipole_estimate(bexp, basis, L, Cov) 0004 % 単一ダイポール推定 0005 % 0006 % Vindx : 推定されたダイポールのインデックス 0007 % : ダイポールの位置座標 = V(Vindx,:) 0008 % J : ダイポール電流強度 0009 % emin : 復元誤差 0010 % 0011 % bexp : MEG-data (Nsensor ,Time) 0012 % basis : Lead Field (Ndipole, Nsensor) 0013 % L = 2 (電流方向) 0014 % Cov : Noise covariance matrix : Cov(sensor,sensor) 0015 % --- 0016 % 0017 % Copyright (C) 2011, ATR All Rights Reserved. 0018 % License : New BSD License(see VBMEG_LICENSE.txt) 0019 0020 if nargin<3, 0021 L = 2; 0022 end; 0023 0024 basis = basis'; % basis(Nsensor, Ndipole) 0025 0026 [N , M] = size(basis); 0027 Npoint = M/L; 0028 0029 if nargin<4, 0030 SS = eye(N,N); 0031 else 0032 SS = pinv(Cov); % Inverce of Noise covariance matrix 0033 end; 0034 0035 err = zeros(Npoint,1); 0036 JJ = zeros(Npoint,L); 0037 0038 id = 0:(L-1); 0039 id = Npoint*id; % [ 0 , Npoint] for L=2 0040 0041 Bnorm = sum(sum(bexp.^2)); 0042 0043 for n = 1:Npoint, 0044 % Lead Field for the n-th dipole 0045 G = basis(: , n + id); 0046 0047 % Min. ( bexp - G*J )^2 0048 GG = inv( G' * SS * G ); % (L x L) matrix inverse 0049 0050 % Diploe current 0051 J = GG * (G' * SS * bexp); 0052 berr = bexp - G * J; 0053 0054 % Reconstruction error 0055 err(n) = sum( berr' * SS * berr )/Bnorm; 0056 JJ(n,:) = sqrt(sum( J.^2, 1)); 0057 end; 0058 0059 [ emin, Vindx ] = min( err ); 0060 J = JJ(Vindx,:);