returns [design] matrix of explanatory variables FORMAT [X Xn IND BF name] = spm_Volterra(SF,BF,name,N); SF{i} - multivariate causes: SF{i}(:,j) = casue i, expansion j BF{i} - Basis functions: BF{i} = basis set for cause i name{i} - name of cause i N - [1 or 2] order of Volterra expansion X - Design Matrix Xn{i} - name of cause i (now including interactions among causes) IND{i} - indices pertaining to cause i (and interactions) ___________________________________________________________________________ For first order expansions spm_Volterra simply convolves the causes (e.g. stick functions) in SF by the basis functions in BF to create a design matrix X. For second order expansions new entries appear in IND, BF and name that correspond to the interaction among the orginal causes (if the events are sufficiently close in time). The basis functions for these are two dimensional and are used to assemble the second order kernel in spm_graph.m. Second order effects are computed for only the first column of SF. ___________________________________________________________________________ @(#)spm_Volterra.m 2.1 Karl Friston 99/05/11
0001 function [X,Xn,IND,BF,name] = spm_Volterra(SF,BF,name,N) 0002 % returns [design] matrix of explanatory variables 0003 % FORMAT [X Xn IND BF name] = spm_Volterra(SF,BF,name,N); 0004 % SF{i} - multivariate causes: SF{i}(:,j) = casue i, expansion j 0005 % BF{i} - Basis functions: BF{i} = basis set for cause i 0006 % name{i} - name of cause i 0007 % N - [1 or 2] order of Volterra expansion 0008 % 0009 % X - Design Matrix 0010 % Xn{i} - name of cause i (now including interactions among causes) 0011 % IND{i} - indices pertaining to cause i (and interactions) 0012 %___________________________________________________________________________ 0013 % 0014 % For first order expansions spm_Volterra simply convolves the causes 0015 % (e.g. stick functions) in SF by the basis functions in BF to create 0016 % a design matrix X. For second order expansions new entries appear 0017 % in IND, BF and name that correspond to the interaction among the 0018 % orginal causes (if the events are sufficiently close in time). 0019 % The basis functions for these are two dimensional and are used to 0020 % assemble the second order kernel in spm_graph.m. Second order effects 0021 % are computed for only the first column of SF. 0022 %___________________________________________________________________________ 0023 % @(#)spm_Volterra.m 2.1 Karl Friston 99/05/11 0024 0025 0026 % Construct X 0027 %=========================================================================== 0028 0029 % 1st order terms 0030 %--------------------------------------------------------------------------- 0031 X = []; 0032 Xn = {}; 0033 IND = cell(1,size(SF,2)); 0034 for i = 1:size(SF,2) 0035 for j = 1:size(BF{i},2) 0036 for k = 1:size(SF{i},2) 0037 x = SF{i}(:,k); 0038 d = 1:length(x); 0039 x = conv(full(x),BF{i}(:,j)); 0040 x = x(d); 0041 X = [X x]; 0042 IND{i} = [IND{i} size(X,2)]; 0043 if size(SF{i},2) > 1 0044 str = [name{i} sprintf('(%i)[%i]',j,(k - 1))]; 0045 else 0046 str = [name{i} sprintf('(%i)',j)]; 0047 end 0048 Xn{end + 1} = str; 0049 end 0050 end 0051 end 0052 0053 % return if first order 0054 %--------------------------------------------------------------------------- 0055 if N == 1, return, end 0056 0057 % 2nd order terms 0058 %--------------------------------------------------------------------------- 0059 k = length(name); 0060 for i = 1:size(SF,2) 0061 for j = i:size(SF,2) 0062 0063 % ensure events can interact 0064 %----------------------------------------------------------- 0065 skip = 0; 0066 if i == j 0067 p = diff(find(SF{i}(:,1))); 0068 skip = (size(BF{i},1) <= min(p)) | ~any(diff(p)); 0069 end 0070 0071 if ~skip 0072 0073 k = k + 1; 0074 ind = []; 0075 bf = {}; 0076 for p = 1:size(BF{i},2) 0077 for q = 1:size(BF{j},2) 0078 ni = [name{i} sprintf('(%i)',p)]; 0079 nj = [name{j} sprintf('(%i)',q)]; 0080 x = SF{i}(:,1); 0081 y = SF{j}(:,1); 0082 x = conv(full(x),BF{i}(:,p)); 0083 y = conv(full(y),BF{j}(:,q)); 0084 x = x(d); 0085 y = y(d); 0086 X = [X x.*y]; 0087 ind = [ind size(X,2)]; 0088 Xn{end + 1} = [ni ' x ' nj]; 0089 bf = [bf {BF{i}(:,p)*BF{j}(:,q)'}]; 0090 end 0091 end 0092 name{k} = [name{i} ' x ' name{j}]; 0093 IND{k} = ind; 0094 BF{k} = bf; 0095 end 0096 end 0097 end