Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Fix issue 395 - pspm_con1 #398

Merged
merged 11 commits into from
Jul 11, 2022
Merged

Fix issue 395 - pspm_con1 #398

Changes from 10 commits
Commits
Show all changes
11 commits
Select commit Hold shift + click to select a range
12479d3
convert nan to 0 in data.stats
Jun 14, 2022
5b04878
Update pspm_con1.m
Jun 14, 2022
9a45ceb
Update pspm_con1.m
Jun 21, 2022
6d589f2
Update pspm_con1.m
Jun 21, 2022
bb9e1aa
Update pspm_con1.m
Jun 21, 2022
72fa735
Update pspm_con1.m
Jun 28, 2022
4fd52a5
Revert "Update pspm_con1.m"
Jun 28, 2022
360202a
Update pspm_con1.m
Jul 5, 2022
4ad7a9e
Update pspm_con1.m
Jul 5, 2022
0c9828a
Update pspm_con1.m
Jul 5, 2022
3ec94e6
Merge branch 'develop' into 395-incorrect-assigned-nan
teddychao Jul 11, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
224 changes: 114 additions & 110 deletions src/pspm_con1.m
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,176 +1,180 @@
function pspm_con1(modelfile, connames, convec, datatype, deletecon, options)
% pspm_con1 creates contrasts on the first level (i.e. within one dataset)
% and saves them to the modelfile to be accessed later
%
% FORMAT:
% pspm_con1 (modelfile, connames, convec, [datatype, deletecon, options])
%
% modelfile: a filename, or cell array of filenames
% connames: a cell array of names for the desired contrasts
% convec: a cell array of contrasts
%
% optional arguments
% datatype: 'param': use all parameter estimates
% 'cond': GLM - contrasts formulated in terms of conditions,
% automatically detects number of basis functions and
% uses only the first one (i.e. without derivatives)
% other models - contrasts based on unique trial names
% 'recon': contrasts formulated in terms of conditions in a GLM,
% reconstructs estimated response from all basis functions
% and uses the peak of the estimated response
% deletecon: should existing contrasts be deleted (1) or
% appended (0)?
% default = 0;
% options: options.zscored: 1 - zscore data
% Restriction: only for non-linear models
% and not when datatype == 'recon'
% 0 - do not zscore data
%__________________________________________________________________________
% PsPM 3.0
% (C) 2008-2015 Dominik R Bach (Wellcome Trust Centre for Neuroimaging)

% $Id$
% $Rev$
% ● Description
% pspm_con1 creates contrasts on the first level (i.e. within one
% dataset) and saves them to the modelfile to be accessed later.
% ● Format
% pspm_con1 (modelfile, connames, convec, [datatype, deletecon, options])
% ● Arguments
% modelfile: a filename, or cell array of filenames.
% connames: a cell array of names for the desired contrasts.
% convec: a cell array of contrasts.
% datatype: an optional structure
% .param: use all parameter estimates
% .cond: GLM - contrasts formulated in terms of conditions,
% automatically detects number of basis functions and uses
% only the first one (i.e. without derivatives) other models
% - contrasts based on unique trial names.
% .recon: contrasts formulated in terms of conditions in a GLM,
% reconstructs estimated response from all basis functions
% and uses the peak of the estimated response
% .deletecon: should existing contrasts be deleted (1) or appended (0)?
% default as 0;
% options:
% .zscored: 1 - zscore data
% Restriction: only for non-linear models
% and not when datatype == 'recon'
% 0 - do not zscore data
% ● Version
% PsPM 3.0
% (C) 2008-2015 Dominik R Bach (Wellcome Trust Centre for Neuroimaging)

%% Initialise
global settings
if isempty(settings)
pspm_init;
end
sts = -1;

% check input arguments
% ------------------------------------------------------------------------
if nargin<1
errmsg=sprintf('No modelfile specified'); warning(errmsg); return;
elseif nargin<2
errmsg=sprintf('No contrast names specified'); warning(errmsg); return;
elseif nargin<3
errmsg=sprintf('No contrasts specified'); warning(errmsg); return;
elseif nargin<4
%% 2 check input arguments
% 2.1 check nargin
if nargin < 1
errmsg = 'No modelfile specified'; warning(errmsg); return;
elseif nargin < 2
errmsg = 'No contrast names specified'; warning(errmsg); return;
elseif nargin < 3
errmsg = 'No contrasts specified'; warning(errmsg); return;
elseif nargin < 4
datatype = 'param';
end;
end
if nargin < 5
deletecon = 0;
end;
end
if nargin < 6
options = struct();
end;

% check & convert filename --
end
% 2.2 check & convert filenames
if ischar(modelfile)
modelfile={modelfile};
modelfile = {modelfile};
elseif ~iscell(modelfile)
warning('Model file must be string or cell array of string.');
end;

% check contrasts --
end
% 2.3 check contrasts
if ~iscell(convec)
warning('Please specify a cell array of contrast vectors.'); return;
end;
warning('Please specify a cell array of contrast vectors.'); return
end
for c=1:numel(convec)
if ~isnumeric(convec{c})
errmsg=sprintf('Contrast #%d is not a numeric vector.', c); warning(errmsg); return;
end;
end;

% check contrast --
if numel(connames)~=numel(convec)
warning('Contrast #%d is not a numeric vector.', c); return
end
end
% 2.4 check contrast
if numel(connames) ~= numel(convec)
warning('Number of contrast names (%d) and number of contrast vectors (%d) don''t match.', ...
numel(connames), numel(convec));
return;
end;

% store for output
return
end
% 2.5 store for output
out_datatype = datatype;
switch datatype
case 'param'
datatype = 'stats';
case {'cond', 'recon'}
otherwise
warning('Unknown datatype');
return;
end;

% set load1_options
return
end
% 2.6 set load1_options
load1_options = struct('zscored',0);
if isfield(options, 'zscored') && options.zscored
load1_options.zscored = 1;
end;

% work on contrasts
% ------------------------------------------------------------------------
for iFn =1:numel(modelfile)
% user output --
end
%% 3 work on contrasts
for iFn = 1:numel(modelfile)
% 3.1 user output --
fprintf('Loading data ... ');

% retrieve stats --
[sts, data, mdltype] = pspm_load1(modelfile{iFn}, datatype, '', load1_options);
if sts == -1
% 3.2 retrieve stats --
[lsts, data, ~] = pspm_load1(modelfile{iFn}, datatype, '', load1_options);
if lsts == -1
warning('ID:invalid_input', 'Could not retrieve stats');
jbrochar marked this conversation as resolved.
Show resolved Hide resolved
return;
end;
% create con structure or retrieve existing contrasts --
end
% 3.3 create con structure or retrieve existing contrasts --
if deletecon == 1
con = []; conno = 0;
else
[sts, con] = pspm_load1(modelfile{iFn}, 'con');
if sts == -1
[lsts, con] = pspm_load1(modelfile{iFn}, 'con');
if lsts == -1
fprintf(' Creating fresh contrast structure.\n');
con = []; conno = 0;
else
conno = numel(con);
end;
end;

% user output --
end
end
% 3.4 user output --
fprintf('\nWriting contrasts to %s\n', modelfile{iFn});
% check number of contrast weights --
% 3.5 check number of contrast weights --
paramno = size(data.stats, 1);
for c=1:numel(convec)
for c = 1:numel(convec)
if numel(convec{c}) > paramno
warning('Contrast (%d) has more weights than statistics (%d) in modelfile %s', ...
jbrochar marked this conversation as resolved.
Show resolved Hide resolved
numel(convec{c}), paramno, modelfile{iFn}); return;
end;
end;
% transform into row vector and right pad with zeros --
conmat=zeros(numel(connames), paramno);
for c=1:numel(convec)
conmat(c,1:numel(convec{c}))=convec{c};
end;

% calculate contrasts --
numel(convec{c}), paramno, modelfile{iFn});
return
end
end
% 3.6 transform into row vector and right pad with zeros --
conmat = zeros(numel(connames), paramno);
for c = 1:numel(convec)
conmat(c,1:numel(convec{c})) = convec{c};
end
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Tip: If all cells of convec are row-vector, I think conmat = cat(1,convec{:}); is similar to your for loop

Copy link
Contributor Author

@teddychao teddychao Jul 5, 2022 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is not written by me. I just checked and I do not think they do have similar features. In the old code convec could be a 2D matrix thus we may need to do for each row..

% 3.7 calculate contrasts
% this automatically replicates contrasts across multiple measures if
% data.stats has more than one column
conval = conmat * data.stats;
% there are issues if data.stats has NaN and the corresponding conmat
% also contains 0
idx_issue = any(isnan(data.stats),2);
l_conmat_r = length(conmat(:,1));
if any(idx_issue(:))
for i_conmat_r = 1:l_conmat_r % if conmat is a 2D matrix, repeat for each row
conmat_array = conmat(i_conmat_r,:);
idx_valid_issue = transpose(conmat_array==0) .* idx_issue;
idx_invalid_issue = transpose(conmat_array~=0) .* idx_issue;
if any(idx_valid_issue) && ~any(idx_invalid_issue)
warning(['Calculated data.stats contain NaNs that are caused by unknown reasons. '...
'data.stats are then used in a contrasts, producing an invalid results (NaN).']);
data_stats_converted = data.stats;
data_stats_converted(isnan(data_stats_converted)) = 0;
conval = conmat * data_stats_converted;
else
warning(['Calculated data.stats contain NaNs that are caused by unknown reasons. '...
'However they were not used in the computation of the contrasts.']);
conval = conmat * data.stats;
end
end
end

% zscored text-output for connames
if isfield(data, 'zscored') && data.zscored
out_zscored = ' (z-scored)';
else
out_zscored = '';
end

% create name matrix if necessary --
% 3.8 create name matrix if necessary --
if size(conval, 2) > 1
for iCon = 1:size(conval, 1)
for iMsr = 1:size(conval, 2)
newconnames{iCon, iMsr} = sprintf('%s - %s%s', connames{iCon}, ...
data.names{iMsr}, out_zscored);
end;
end;
end
end
else
newconnames = connames;
end;

% save contrasts
for iCon=1:numel(conval)
end
% 3.9 save contrasts
for iCon = 1:numel(conval)
con(conno+iCon).type = out_datatype;
con(conno+iCon).name = newconnames{iCon};
con(conno+iCon).con = conval(iCon);
indx = mod(iCon, size(conval, 1));
if indx == 0, indx = size(conval, 1); end;
if indx == 0, indx = size(conval, 1); end
con(conno+iCon).convec = conmat(indx, :);
end;
end
pspm_load1(modelfile{iFn}, 'savecon', con);
end;
end