[DOC] improve stats doc and warnings

.github/dependabot.yml

@@ -6,5 +6,3 @@ updates:
     directory: "/"
     schedule:
       interval: "weekly"
-    labels:
-      - "exclude-from-changelog"

demos/MoAE/moae_01_bids_app.m

@@ -48,7 +48,7 @@
 cpp_spm();
 
 %% Dopwnload the dataset
-download_data = false;
+download_data = true;
 clean = false;
 download_moae_ds(download_data, clean);
 

docs/source/FAQ.md

@@ -6,12 +6,12 @@
 
 ### How can I prevent from having SPM windows pop up all the freaking time?
 
-Running large number of batches when the GUI of Matlab is active can be
+Running large number of batches when the GUI of MATLAB is active can be
 annoying, as SPM windows will always pop up and become active instead of running
 in the background like most users would prefer to.
 
 One easy solution is to add a `spm_my_defaults` function with the following
-content in the Matlab path, or in the directory where you are running your
+content in the MATLAB path, or in the directory where you are running your
 scripts or command from.
 
 ```matlab
@@ -30,9 +30,9 @@ SPM runs in command line mode.
 ### How can I run any of this from the command line?
 
 Related to the previous question but more radical. You can run most analysis
-from within your terminal without starting the Matlab graphic interface.
+from within your terminal without starting the MATLAB graphic interface.
 
-For this you first need to know where is the Matlab application. Here are the
+For this you first need to know where is the MATLAB application. Here are the
 typical location depending on your operating system (where `XXx` corresponds to
 the version you use).
 
@@ -183,7 +183,7 @@ it.
 
 Ideally you would want to add a JSON file to add metadata about those ROIs.
 
-You can use bids matlab to help you create bids valid filenames.
+You can use bids-matlab to help you create BIDS valid filenames.
 
 ```matlab
   >> name_spec.ext = '.nii';
@@ -213,7 +213,7 @@ This is similar to the way you can "select" only certain files to preprocess
 with [fmriprep](https://fmriprep.org/en/stable/faq.html#how-do-i-select-only-certain-files-to-be-input-to-fmriprep).
 
 You can use a `opt.bidsFilterFile` field in your options
-to define a typical images "bold", "T1w" in terms of their bids entities.
+to define a typical images "bold", "T1w" in terms of their BIDS entities.
 The default value is:
 
 ```matlab
@@ -256,10 +256,10 @@ you would define this file like this:
 
 You can select a subset of your data by using the `opt.query`.
 
-This will create a "filter" that bids matlab will use to only "query" and
+This will create a "filter" that bids-matlab will use to only "query" and
 retrieve the subset of files that match the requirement of that filter
 
-In "pure" bids matlab it would look like:
+In "pure" bids-matlab it would look like:
 
 ```matlab
   BIDS = bids.layout(path_to_my_dataset)
@@ -345,7 +345,7 @@ If your data was preprocessed with fmriprep, cpp_spm will first need to copy,
 unzip and smooth the data into a `cpp_spm-preproc` folder
 
 Here is an example of how the data is organized for the MoAE fmriprep demo and
-what the `cpp_spm` bids call would look like.
+what the `cpp_spm` BIDS call would look like.
 
 ```bash
 ├── inputs

src/IO/getData.m

@@ -44,15 +44,14 @@
     printToScreen(msg, opt);
   end
 
-  msg = sprintf('Getting data from:\n %s\n', pathToPrint(bidsDir));
-  printToScreen(msg, opt);
-
   validationInputFile(bidsDir, 'dataset_description.json');
 
-  BIDS = bids.layout(bidsDir, 'use_schema', opt.useBidsSchema);
+  BIDS = bids.layout(bidsDir, 'use_schema', opt.useBidsSchema, 'verbose', opt.verbosity > 0);
 
   if strcmp(opt.pipeline.type, 'stats')
     if exist(fullfile(opt.dir.raw, 'layout.mat'), 'file') == 2
+      msg = sprintf('\nLoading BIDS raw layout from:\n\t%s\n', fullfile(opt.dir.raw, 'layout.mat'));
+      printToScreen(msg, opt);
       tmp = load(fullfile(opt.dir.raw, 'layout.mat'), 'BIDS');
       if isempty(fieldnames(tmp))
         BIDS.raw = bids.layout(opt.dir.raw);

src/IO/overwriteDir.m

@@ -7,7 +7,7 @@ function overwriteDir(directory, opt)
   % (C) Copyright 2021 CPP_SPM developers
 
   if exist(directory, 'dir') == 7
-    msg = sprintf('overwriting directory: %s \n', pathToPrint(directory));
+    msg = sprintf('\noverwriting directory:\n\t%s\n\n', pathToPrint(directory));
     errorHandling(mfilename(), 'overWritingDir', msg, true, opt.verbosity);
     rmdir(directory, 's');
   end

src/batches/stats/setBatchSubjectLevelGLMSpec.m

@@ -52,24 +52,26 @@
   fmri_spec.timing.units = 'secs';
   fmri_spec.timing.RT = TR;
 
+  % unique is used in case data was acquired multiband
   nbTimeBins = numel(unique(sliceOrder));
   fmri_spec.timing.fmri_t = nbTimeBins;
 
   % If no reference slice is given for STC,
   % then STC took the mid-volume as reference time point for the GLM.
   % When no STC was done, this is usually a good way to do it too.
   if isempty(opt.stc.referenceSlice)
-    refBin = floor(nbTimeBins / 2);
+    refBin = nbTimeBins / 2;
   else
     refBin = opt.stc.referenceSlice / TR;
   end
+  refBin = floor(refBin);
   fmri_spec.timing.fmri_t0 = refBin;
 
   % Create ffxDir if it doesnt exist
   % If it exists, issue a warning that it has been overwritten
   ffxDir = getFFXdir(subLabel, opt);
   overwriteDir(ffxDir, opt);
-  printToScreen(sprintf(' output dir: %s\n', pathToPrint(ffxDir)), opt);
+  printToScreen(sprintf('\n output dir:\n\t%s\n\n', pathToPrint(ffxDir)), opt);
   fmri_spec.dir = {ffxDir};
 
   fmri_spec.fact = struct('name', {}, 'levels', {});
@@ -177,25 +179,52 @@
 
   filter = fileFilterForBold(opt, subLabel);
 
-  if ~opt.stc.skip
+  [name, version] = generatedBy(BIDS);
+  tokens = strsplit(version);
+  % TODO implement differently for fmriprep >=20.2.4
+  % https://fmriprep.org/en/stable/changes.html#october-04-2021
+  skip = opt.stc.skip || ...
+        (strcmp(name, 'fMRIPrep') && ...
+         (str2num(tokens{1}) > 20 || ...
+          str2num(tokens{1}) == 20 && str2num(tokens{2}) >= 4));  %#ok<*ST2NM>
+
+  sliceOrder = [];
+  if ~skip
     % Get slice timing information.
     % Necessary to make sure that the reference slice used for slice time
     % correction is the one we center our model on;
+
+    % temporary silence warnings and only throw a single warning
+    % after that if necessary
+    oldVerbosity = opt.verbosity;
+    opt.verbosity = 0;
+
     sliceOrder = getAndCheckSliceOrder(BIDS, opt, filter);
-  else
-    sliceOrder = [];
+
+    opt.verbosity = oldVerbosity;
   end
 
   if isempty(sliceOrder) && ~opt.dryRun
-    % no slice order defined here (or different across tasks)
-    % so we fall back on using the number of
-    % slice in the first bold image to set the number of time bins
-    % we will use to upsample our model during regression creation
+
     fileName = bids.query(BIDS, 'data', filter);
     hdr = spm_vol(fileName{1});
 
-    % we are assuming axial acquisition here
+    % TODO we are assuming axial acquisition here
     sliceOrder = 1:hdr(1).dim(3);
+
+    wng = ['\n\n', ...
+           'Slice timing information was missing for at least one run,\n', ...
+           'or was inconsistent across runs.', ...
+           '\n', ...
+           'Will be using the number of slices as the number of bins\n', ...
+           'for temporal upsampling before convolution.', ...
+           '\n', ...
+           'If your data was processed with fMRIprep < 20.2.4, this is expected.\n'];
+    % note that with multiband
+    % this may lead to more time bins that used in reality at acquisition
+
+    errorHandling(mfilename(), 'noSliceTimingInfoForGlm', wng, true, opt.verbosity);
+
   end
 
 end

src/bids_model/BidsModel.m

@@ -106,7 +106,10 @@
           isempty(model.Options.HighPassFilterCutoffHz)
 
         obj.bidsModelError('noHighPassFilter', ...
-                           sprintf('No high-pass filter for node "%s"', nodeName));
+                           sprintf(['No high-pass filter for node "%s"\n', ...
+                                    'Will use default falue instead: %f.'], ...
+                                   nodeName, ...
+                                   spm_get_defaults('stats.fmri.hpf')));
 
       else
 
@@ -147,7 +150,8 @@
         HRFderivatives = model.HRF.Model;
       catch
         obj.bidsModelError('noHRFderivatives', ...
-                           sprintf('No HRF derivatives for node "%s"', nodeName));
+                           sprintf(['No HRF derivatives for node "%s".\n', ...
+                                    'Will use SPM canonical HRF.\n'], nodeName));
       end
 
       HRFderivatives = lower(strrep(HRFderivatives, ' ', ''));
@@ -187,7 +191,7 @@
         mask = model.Options.Mask;
       catch
         msg = sprintf(['No mask for Node "%s".', ...
-                       '\nSpecify one in "Node.Options.Mask"'], ...
+                       '\nSpecify one in "Node.Options.Mask".\n'], ...
                       nodeName);
         obj.bidsModelError('noMask', msg);
       end
@@ -213,7 +217,10 @@
       else
 
         obj.bidsModelError('noInclMaskThresh', ...
-                           sprintf('No inclusive mask threshold for Node "%s"', nodeName));
+                           sprintf(['No inclusive mask threshold for Node "%s".\n', ...
+                                    'Will use default falue instead: %f.'], ...
+                                   nodeName, ...
+                                   spm_get_defaults('mask.thresh')));
 
       end
 
@@ -241,13 +248,17 @@
 
       else
         obj.bidsModelError('noTemporalCorrection', ...
-                           sprintf('No temporal correlation correction for Node "%s"', nodeName));
+                           sprintf(['No temporal correlation correction for Node "%s"\n', ...
+                                    'Will use default falue instead: %s.'], ...
+                                   nodeName, ...
+                                   spm_get_defaults('stats.fmri.cvi')));
+
       end
 
     end
 
     function bidsModelError(obj, id, msg)
-      msg = sprintf('\n\n%s in BIDS stats model named "%s"\n', msg, obj.Name);
+      msg = sprintf('\n\nFor BIDS stats model named: "%s"\n%s\n', obj.Name, msg);
       errorHandling(mfilename(), id, msg, obj.tolerant, obj.verbose);
     end
 

src/defaults/checkOptions.m

@@ -108,7 +108,8 @@
   %
   %
   %     - ``opt.stc.skip = false`` - boolean flag to skip slice time correction or not.
-  %     - ``opt.stc.referenceSlice = []`` - reference slice for the slice timing correction.
+  %     - ``opt.stc.referenceSlice = []`` - reference slice (in seconds)
+  %       for the slice timing correction.
   %       If left emtpy the mid-volume acquisition time point will be selected at run time.
   %     - ``opt.stc.sliceOrder = []`` - To be used if SPM can't extract slice info. NOT RECOMMENDED,
   %       if you know the order in which slices were acquired, you should be able to recompute

src/messages/printProcessingSubject.m

@@ -7,8 +7,7 @@ function printProcessingSubject(iSub, subLabel, opt)
   %
   % (C) Copyright 2019 CPP_SPM developers
 
-  msg = sprintf([ ...
-                 ' PROCESSING SUBJECT No.: %i ' ...
+  msg = sprintf([' PROCESSING SUBJECT No.: %i ' ...
                  'SUBJECT LABEL : %s \n'], ...
                 iSub, subLabel);
 

src/workflows/stats/bidsFFX.m

@@ -180,6 +180,9 @@ function checkRootNode(opt)
     errorHandling(mfilename(), 'tooManySpaces', msg, false, opt.verbosity);
   end
 
+  msg = sprintf('\n PROCESSING NODE: %s\n', thisNode.Name);
+  printToScreen(msg, opt, 'format', '*blue');
+
 end
 
 function status = subjectHasData(BIDS, opt, subLabel)

src/workflows/stats/bidsRFX.m

@@ -103,6 +103,9 @@
 
       for i = 1:numel(datasetNodes)
 
+        msg = sprintf('\n PROCESSING NODE: %s\n', nodeName);
+        printToScreen(msg, opt, 'format', '*blue');
+
         matlabbatch = {};
 
         nodeName = datasetNodes{i}.Name;

src/workflows/stats/bidsResults.m

@@ -186,6 +186,9 @@
       continue
     end
 
+    msg = sprintf('\n PROCESSING NODE: %s\n', node.Name);
+    printToScreen(msg, opt, 'format', '*blue');
+
     % Depending on the level step we migh have to define a matlabbatch
     % for each subject or just on for the whole group
     switch lower(node.Level)