-
Notifications
You must be signed in to change notification settings - Fork 29
/
mTRFtransform.m
182 lines (157 loc) · 6.11 KB
/
mTRFtransform.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
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
function fmodel = mTRFtransform(bmodel,resp,varargin)
%MTRFTRANSFORM Transform a backward model into a forward model.
% FMODEL = MTRFTRANSFORM(BMODEL,RESP) transforms a backward decoding
% model (stimulus to neural response) into the corresponding forward
% encoding model (neural response to stimulus) as per Haufe et al.
% (2014), enabling the neurophysiological interpretation of the backward
% model weights. RESP are the neural responses that were used to compute
% the original backward model.
%
% MTRFTRANSFORM returns the transformed model in a structure with the
% following fields:
% 'w' -- transformed model weights (xvar-by-nlag-by-yvar)
% 't' -- time lags (ms)
% 'fs' -- sample rate (Hz)
% 'Dir' -- direction of causality (forward=1, backward=-1)
% 'type' -- type of model (multi-lag, single-lag)
%
% If RESP is a matrix, it is assumed that the rows correspond to
% observations and the columns to variables, unless otherwise stated via
% the 'dim' parameter (see below). If RESP is a vector, it is assumed
% that the first non-singleton dimension corresponds to observations.
%
% If RESP is a cell array containing multiple trials, the covariance
% matrices of each trial are summed to produce one transformed model.
%
% [...] = MTRFTRANSFORM(...,'PARAM1',VAL1,'PARAM2',VAL2,...) specifies
% additional parameters and their values. Valid parameters are the
% following:
%
% Parameter Value
% 'dim' A scalar specifying the dimension to work along: pass
% in 1 to work along the columns (default) or 2 to work
% along the rows.
% 'split' A scalar specifying the number of segments in which to
% split each trial of data when computing the covariance
% matrices. This is useful for reducing memory usage on
% large datasets. By default, the entire trial is used.
% 'zeropad' A numeric or logical specifying whether to zero-pad the
% outer rows of the design matrix or delete them: pass in
% 1 to zero-pad them (default), or 0 to delete them.
% 'verbose' A numeric or logical specifying whether to execute in
% verbose mode: pass in 1 for verbose mode (default), or
% 0 for non-verbose mode.
%
% See also MTRFTRAIN, MTRFMULTITRAIN.
%
% mTRF-Toolbox https://github.com/mickcrosse/mTRF-Toolbox
% References:
% [1] Haufe S, Meinecke F, Gorgen K, Dahne S, Haynes JD, Blankertz B,
% Bießmann F (2014) On the interpretation of weight vectors of
% linear models in multivariate neuroimaging. NeuroImage
% 87:96-110.
% [2] Crosse MC, Di Liberto GM, Bednar A, Lalor EC (2016) The
% multivariate temporal response function (mTRF) toolbox: a MATLAB
% toolbox for relating neural signals to continuous stimuli. Front
% Hum Neurosci 10:604.
% Authors: Mick Crosse <[email protected]>
% Adam Bednar <[email protected]>
% Emily Teoh <[email protected]>
% Giovanni Di Liberto <[email protected]>
% Copyright 2014-2024 Lalor Lab, Trinity College Dublin.
% Parse input arguments
arg = parsevarargin(varargin);
% Format data in cell arrays
[resp,nobs] = formatcells(resp,arg.dim,arg.split);
nfold = numel(resp);
% Convert time lags to samples
tmin = bmodel.t(1)*bmodel.fs/1e3;
tmax = bmodel.t(end)*bmodel.fs/1e3;
% Predict output
pred = mTRFpredict([],resp,bmodel,'dim',arg.dim,'zeropad',arg.zeropad,...
'verbose',0);
% Convert to cell array
if ~iscell(pred)
pred = {pred};
end
% Truncate output
if ~arg.zeropad
resp = truncate(resp,tmin,tmax,nobs);
end
% Verbose mode
if arg.verbose
v = verbosemode([],0,nfold);
end
% Initialize variables
Cxx = 0;
Css = 0;
for i = 1:nfold
% Compute covariance matrices
Cxx = Cxx + resp{i}'*resp{i};
Css = Css + pred{i}'*pred{i};
% Verbose mode
if arg.verbose
v = verbosemode(v,i,nfold);
end
end
% Transform backward model weights
for i = 1:length(Css)
bmodel.w(:,:,i) = Cxx*bmodel.w(:,:,i)/Css(i,i);
end
bmodel.w = permute(bmodel.w,[3,2,1]);
% Format output arguments
fmodel = struct('w',fliplr(bmodel.w),'t',-fliplr(bmodel.t),...
'fs',bmodel.fs,'Dir',-bmodel.Dir,'type',bmodel.type);
% Verbose mode
if arg.verbose
verbosemode(v,nfold+1,nfold,fmodel);
end
function v = verbosemode(v,fold,nfold,model)
%VERBOSEMODE Execute verbose mode.
% V = VERBOSEMODE(V,FOLD,NFOLD,MODEL) prints details about the progress
% of the main function to the screen.
if fold == 0
v = struct('msg',[],'h',[],'tocs',0);
fprintf('\nTransform decoding model\n')
fprintf('Computing covariance matrices\n')
v.msg = ['%d/%d [',repmat(' ',1,nfold),']\n'];
v.h = fprintf(v.msg,fold,nfold);
elseif fold <= nfold
if fold == 1 && toc < 0.1
pause(0.1)
end
v.tocs = v.tocs + toc;
fprintf(repmat('\b',1,v.h))
v.msg = ['%d/%d [',repmat('=',1,fold),repmat(' ',1,nfold-fold),'] - ',...
'%.3fs/fold\n'];
v.h = fprintf(v.msg,fold,nfold,v.tocs/fold);
tic
end
if fold == nfold
fprintf('Transforming model');
elseif fold > nfold
fprintf(' - %.3fs\n',toc)
modelsummary(model)
end
function arg = parsevarargin(varargin)
%PARSEVARARGIN Parse input arguments.
% [PARAM1,PARAM2,...] = PARSEVARARGIN('PARAM1',VAL1,'PARAM2',VAL2,...)
% parses the input arguments of the main function.
% Create parser object
p = inputParser;
% Dimension to work along
errorMsg = 'It must be a positive integer scalar within indexing range.';
validFcn = @(x) assert(x==1||x==2,errorMsg);
addParameter(p,'dim',1,validFcn);
% Split data
errorMsg = 'It must be a positive integer scalar.';
validFcn = @(x) assert(isnumeric(x)&&isscalar(x),errorMsg);
addParameter(p,'split',1,validFcn);
% Boolean arguments
errorMsg = 'It must be a numeric scalar (0,1) or logical.';
validFcn = @(x) assert(x==0||x==1||islogical(x),errorMsg);
addParameter(p,'zeropad',true,validFcn); % zero-pad design matrix
addParameter(p,'verbose',true,validFcn); % verbose mode
% Parse input arguments
parse(p,varargin{1,1}{:});
arg = p.Results;