Compute an example of input signals to an MAR model
(input computed by subtracting current sources with one-step prediction
from the estimated current sources)
Input
MAR: MAR matrix [Nvact x Nvact, sparse double]
Delta: Time lag matrix [Nvact x Nvact, sparse double]
Z: Estimated current sources [1 x 1, cell][Nvact x Nt, double]
Output
omega: Input signals to an MAR model [Nvact x Nt, double]
2015/10/07 M.Fukushima
Copyright (C) 2011, ATR All Rights Reserved.
License : New BSD License(see VBMEG_LICENSE.txt)0001 function [omega] = calc_MARinput(MAR, Delta, Z) 0002 % Compute an example of input signals to an MAR model 0003 % (input computed by subtracting current sources with one-step prediction 0004 % from the estimated current sources) 0005 % 0006 % Input 0007 % MAR: MAR matrix [Nvact x Nvact, sparse double] 0008 % Delta: Time lag matrix [Nvact x Nvact, sparse double] 0009 % Z: Estimated current sources [1 x 1, cell][Nvact x Nt, double] 0010 % 0011 % Output 0012 % omega: Input signals to an MAR model [Nvact x Nt, double] 0013 % 0014 % 2015/10/07 M.Fukushima 0015 % 0016 % Copyright (C) 2011, ATR All Rights Reserved. 0017 % License : New BSD License(see VBMEG_LICENSE.txt) 0018 0019 % Estimated current sources 0020 Jfilt = Z{1}; 0021 0022 % # of sources 0023 Nvact = size(MAR,1); 0024 % # of time instances 0025 Nt = size(Jfilt,2); 0026 % Binarized structural connectivity matrix 0027 Ind = (Delta~=0); 0028 0029 % Connection numbers 0030 sum_ind = sum(Ind,2); 0031 % Time lags 0032 unidlt = full(setdiff(unique(Delta),0)); 0033 0034 % Store matrix entries into a smaller size of matrix 0035 Indx = zeros(Nvact,max(sum_ind)); 0036 Deltax = zeros(Nvact,max(sum_ind)); 0037 MARx = zeros(Nvact,max(sum_ind)); 0038 for nv = 1:Nvact 0039 Indx(nv,1:sum_ind(nv)) = find(Ind(nv,:)==1); 0040 Deltax(nv,1:sum_ind(nv)) = Delta(nv,Delta(nv,:)~=0); 0041 MARx(nv,1:sum_ind(nv)) = MAR(nv,Ind(nv,:)~=0); 0042 end 0043 0044 % One step forward prediction with an MAR model 0045 Jpred = zeros(Nvact, Nt); 0046 for t = 1:Nt 0047 Jtmp = [zeros(Nvact, unidlt(end)) Jfilt]; 0048 t_end = t + unidlt(end); 0049 0050 % Prediction [MEX ver.] 0051 AJpred = atimesj(MARx, Indx, Deltax, sum_ind, ... 0052 Jtmp(:,t_end:-1:(t_end-unidlt(end)))); 0053 % % Prediction [Matlab ver.] 0054 % AJpred = zeros(Nvact,1); 0055 % for nv = 1:Nvact 0056 % for ii = 1:sum_ind(nv) 0057 % AJpred(nv) = AJpred(nv) + ... 0058 % MARx(nv,ii) * Jtmp(Indx(nv,ii),t_end-Deltax(nv,ii)); 0059 % end 0060 % end 0061 0062 % Current sources of one step prediction 0063 Jpred(:,t) = AJpred; 0064 end 0065 0066 % Input signals to an MAR model 0067 omega = Jfilt - Jpred;