From bb2aab1b0ff037bc48bec56c39d215d34f5bfeb7 Mon Sep 17 00:00:00 2001 From: Remi Gau Date: Wed, 31 Jul 2024 12:16:34 +0200 Subject: [PATCH 1/7] start adding QA workflow --- .github/workflows/tests.yml | 5 + demos/openneuro/Makefile | 30 ++-- demos/openneuro/ds000114_preproc_run.m | 27 ++-- docs/source/API.rst | 2 +- docs/source/QA.rst | 2 +- src/IO/getData.m | 21 +-- src/utils/setUpWorkflow.m | 10 +- src/workflows/QA/bidsQA.m | 37 +++++ .../{bidsQApreproc.m => QA/bidsQAbidspm.m} | 6 +- .../mriqcQA.m => workflows/QA/bidsQAmriqc.m} | 0 src/workflows/QA/bidsQApreproc.m | 129 ++++++++++++++++++ src/workflows/preproc/bidsResliceTpmToFunc.m | 2 +- src/workflows/preproc/bidsSpatialPrepro.m | 2 +- tests/tests_workflows/QA/test_bidsQA.m | 80 +++++++++++ .../QA/test_bidsQAmriqc.m} | 4 +- 15 files changed, 318 insertions(+), 39 deletions(-) create mode 100644 src/workflows/QA/bidsQA.m rename src/workflows/{bidsQApreproc.m => QA/bidsQAbidspm.m} (98%) rename src/{QA/mriqcQA.m => workflows/QA/bidsQAmriqc.m} (100%) create mode 100644 src/workflows/QA/bidsQApreproc.m create mode 100644 tests/tests_workflows/QA/test_bidsQA.m rename tests/{tests_QA/test_mriqcQA.m => tests_workflows/QA/test_bidsQAmriqc.m} (81%) diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index 3b5f20dd0..818bac6c6 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -117,6 +117,11 @@ jobs: unzip download mv moae_fmriprep fmriprep + - name: Get data for testing QA + run: | + cd demos/openneuro/ + make data_ds000114_mriqc + - name: Prepare test data unix run: | cd tests diff --git a/demos/openneuro/Makefile b/demos/openneuro/Makefile index f714be174..b540e3c55 100644 --- a/demos/openneuro/Makefile +++ b/demos/openneuro/Makefile @@ -14,13 +14,15 @@ clean: data: data_ds000001 data_ds000114 data_ds001168 -data_ds000001: +data_ds000001_fmriprep: + cd inputs && datalad install ///openneuro-derivatives/ds000001-fmriprep + cd inputs/ds000001-fmriprep && datalad get sub-*/func/*tsv -J 12 + cd inputs/ds000001-fmriprep && datalad get sub-*/func/*json -J 12 + +data_ds000001: data_ds000001_fmriprep mkdir -p inputs cd inputs && datalad install ///openneuro/ds000001 - cd inputs && datalad install ///openneuro-derivatives/ds000001-fmriprep cd inputs/ds000001 && datalad get sub-0[1-5] -J 3 - cd inputs/ds000001-fmriprep && datalad get sub-0[1-5]/func/*tsv -J 12 - cd inputs/ds000001-fmriprep && datalad get sub-0[1-5]/func/*json -J 12 cd inputs/ds000001-fmriprep && datalad get sub-0[1-5]/anat/*MNI*desc-preproc*.nii.gz -J 12 cd inputs/ds000001-fmriprep && datalad get sub-0[1-5]/func/*MNI*desc-preproc*.nii.gz -J 12 cd inputs/ds000001-fmriprep && datalad get sub-0[1-3]/func/*MNI*desc-*bold.nii.gz -J 12 @@ -28,16 +30,24 @@ data_ds000001: data_ds000114_raw: mkdir -p inputs cd inputs && datalad install ///openneuro/ds000114 - cd inputs/ds000114 && datalad get sub-0[1-2]/ses-*/anat/*T1w*.nii.gz -J 12 - cd inputs/ds000114 && datalad get sub-0[1-2]/ses-*/func/*linebisection* -J 12 + cd inputs/ds000114 && datalad get sub-0[1-5]/ses-*/anat/*T1w*.nii.gz -J 12 + cd inputs/ds000114 && datalad get sub-0[1-5]/ses-*/func/*bold* -J 12 -data_ds000114: +data_ds000114_mriqc: + mkdir -p inputs + cd inputs && datalad install ///openneuro-derivatives/ds000114-mriqc + cd inputs/ds000114-mriqc && datalad get *.tsv + +data_ds000114_fmriprep: mkdir -p inputs - cd inputs && datalad install ///openneuro/ds000114 cd inputs && datalad install ///openneuro-derivatives/ds000114-fmriprep + cd inputs/ds000114-fmriprep && datalad get sub-0*/ses-*/func/*tsv -J 12 + cd inputs/ds000114-fmriprep && datalad get sub-0*/ses-*/func/*json -J 12 + +data_ds000114: data_ds000114_fmriprep + mkdir -p inputs + cd inputs && datalad install ///openneuro/ds000114 cd inputs/ds000114-fmriprep && datalad get sub-0[1-2]/anat/*MNI*desc-preproc*.nii.gz -J 12 - cd inputs/ds000114-fmriprep && datalad get sub-0[1-2]/ses-*/func/*tsv -J 12 - cd inputs/ds000114-fmriprep && datalad get sub-0[1-2]/ses-*/func/*json -J 12 cd inputs/ds000114-fmriprep && datalad get sub-0[1-2]/ses-*/func/*MNI*_mask.nii.gz -J 12 cd inputs/ds000114-fmriprep && datalad get sub-0[1-2]/ses-*/func/*MNI*desc-preproc*bold.nii.gz -J 12 diff --git a/demos/openneuro/ds000114_preproc_run.m b/demos/openneuro/ds000114_preproc_run.m index 72e9bb91e..7e6be75c4 100644 --- a/demos/openneuro/ds000114_preproc_run.m +++ b/demos/openneuro/ds000114_preproc_run.m @@ -1,3 +1,5 @@ +% Show how to use filter to only include one session to preprocess. + % (C) Copyright 2023 bidspm developers clear; @@ -6,24 +8,29 @@ addpath(fullfile(pwd, '..', '..')); bidspm(); -participant_label = {'01'}; +participant_label = {'04', '05'}; TASK = 'linebisection'; +session_to_select = {'test', 'retest'}; root_dir = fileparts(mfilename('fullpath')); bids_dir = fullfile(root_dir, 'inputs', 'ds000114'); output_dir = fullfile(root_dir, 'outputs', 'ds000114', 'derivatives'); preproc = fullfile(output_dir, 'bidspm-preproc'); -bidspm(bids_dir, output_dir, 'subject', ... - 'participant_label', participant_label, ... - 'action', 'copy', ... - 'task', TASK, ... - 'skip_validation', true, ... - 'verbosity', 3); +%% Copy +% To make sure we still have copied all the data +% bidspm(bids_dir, output_dir, 'subject', ... +% 'participant_label', participant_label, ... +% 'action', 'copy', ... +% 'skip_validation', true, ... +% 'verbosity', 3); -bids_filter_file = struct( ... - 'bold', struct('modality', 'func', 'suffix', 'bold', 'ses', 'retest'), ... - 't1w', struct('modality', 'anat', 'suffix', 'T1w')); +%% Set filter and preprocess +bids_filter_file = struct('bold', struct('modality', 'func', ... + 'suffix', 'bold', ... + 'ses', {session_to_select}), ... + 't1w', struct('modality', 'anat', ... + 'suffix', 'T1w')); bidspm(bids_dir, output_dir, 'subject', ... 'participant_label', participant_label, ... 'action', 'preprocess', ... diff --git a/docs/source/API.rst b/docs/source/API.rst index f03955971..4f16e8e37 100644 --- a/docs/source/API.rst +++ b/docs/source/API.rst @@ -33,7 +33,7 @@ workflows .. autofunction:: src.workflows.bidsCheckVoxelSize .. autofunction:: src.workflows.bidsCopyInputFolder .. autofunction:: src.workflows.bidsInverseNormalize -.. autofunction:: src.workflows.bidsQApreproc +.. autofunction:: src.workflows.bidsQAbidspm .. autofunction:: src.workflows.bidsRename .. autofunction:: src.workflows.bidsReport diff --git a/docs/source/QA.rst b/docs/source/QA.rst index da704f1ca..6cb0b016c 100644 --- a/docs/source/QA.rst +++ b/docs/source/QA.rst @@ -10,7 +10,7 @@ Quality control - :func:`anatQA` -- :func:`bidsQApreproc` +- :func:`bidsQAbidspm` .. _fig_spatialPrepro-reports: .. figure:: preprocessing/images/bidsSpatialPrepro/out_report.png diff --git a/src/IO/getData.m b/src/IO/getData.m index a770de80b..d49f90b98 100644 --- a/src/IO/getData.m +++ b/src/IO/getData.m @@ -8,19 +8,20 @@ % % :param opt: Options chosen for the analysis. % See :func:`checkOptions`. - % % :type opt: structure % % :param bidsDir: the directory where the data is ; defaults to % ``fullfile(opt.dataDir, '..', 'derivatives', 'bidspm')`` - % % :type bidsDir: char % % :param indexDependencies: Use ``'index_dependencies', true`` % in bids.layout. - % % :type indexDependencies: logical % + % :param layoutFilter: filter to pass to bids.layout. default = struct() + % :type layoutFilter: struct + % + % % :return: opt % :rtype: structure % @@ -38,6 +39,7 @@ addRequired(args, 'opt', @isstruct); addRequired(args, 'bidsDir', isFolder); addOptional(args, 'indexDependencies', true, @islogical); + addOptional(args, 'layoutFilter', struct(), @isstruct); try parse(args, varargin{:}); @@ -60,6 +62,8 @@ indexDependencies = false; end + layoutFilter = args.Results.layoutFilter; + anatOnly = false; if isfield(opt, 'anatOnly') anatOnly = opt.anatOnly; @@ -72,19 +76,18 @@ validationInputFile(bidsDir, 'dataset_description.json'); - layout_filter = struct([]); if ~isempty(opt.subjects{1}) && ~ismember('', opt.subjects) - layout_filter = struct('sub', {opt.subjects}); + layoutFilter.sub = opt.subjects; end if anatOnly - layout_filter(1).modality = {'anat'}; + layoutFilter(1).modality = {'anat'}; end BIDS = bids.layout(bidsDir, ... 'use_schema', opt.useBidsSchema, ... 'verbose', opt.verbosity > 1, ... - 'filter', layout_filter, ... + 'filter', layoutFilter, ... 'index_dependencies', indexDependencies); if strcmp(opt.pipeline.type, 'stats') && ~opt.pipeline.isBms @@ -96,7 +99,7 @@ if isempty(fieldnames(tmp)) BIDS.raw = bids.layout(opt.dir.raw, ... - 'filter', layout_filter, ... + 'filter', layoutFilter, ... 'index_dependencies', indexDependencies); else BIDS.raw = tmp.BIDS; @@ -104,7 +107,7 @@ else BIDS.raw = bids.layout(opt.dir.raw, ... 'verbose', opt.verbosity > 1, ... - 'filter', layout_filter); + 'filter', layoutFilter); end end diff --git a/src/utils/setUpWorkflow.m b/src/utils/setUpWorkflow.m index c7aa3aaf5..f57745357 100644 --- a/src/utils/setUpWorkflow.m +++ b/src/utils/setUpWorkflow.m @@ -9,7 +9,7 @@ % % USAGE:: % - % [BIDS, opt, group] = setUpWorkflow(opt, workflowName, bidsDir, indexData) + % [BIDS, opt] = setUpWorkflow(opt, workflowName, bidsDir, indexData, index_dependencies) % % :param opt: Options chosen for the analysis. % See :func:`checkOptions`. @@ -32,6 +32,10 @@ % in bids.layout. default = true % :type index_dependencies: logical % + % :param layoutFilter: filter to pass to bids.layout. default = struct() + % :type layoutFilter: struct + % + % % :returns: ``BIDS`` layout returned by ``getData``, ``opt`` options checked % :rtype: structure, structure % @@ -47,6 +51,7 @@ addOptional(args, 'bidsDir', '', @ischar); addOptional(args, 'indexData', true, @islogical); addOptional(args, 'indexDependencies', true, @islogical); + addOptional(args, 'layoutFilter', struct(), @isstruct); parse(args, varargin{:}); @@ -55,6 +60,7 @@ bidsDir = args.Results.bidsDir; indexData = args.Results.indexData; indexDependencies = args.Results.indexDependencies; + layoutFilter = args.Results.layoutFilter; if isempty(bidsDir) bidsDir = opt.dir.input; @@ -67,7 +73,7 @@ opt = loadAndCheckOptions(opt); if indexData - [BIDS, opt] = getData(opt, bidsDir, indexDependencies); + [BIDS, opt] = getData(opt, bidsDir, indexDependencies, layoutFilter); end if strcmp(opt.pipeline.type, 'stats') && ... diff --git a/src/workflows/QA/bidsQA.m b/src/workflows/QA/bidsQA.m new file mode 100644 index 000000000..13e283b35 --- /dev/null +++ b/src/workflows/QA/bidsQA.m @@ -0,0 +1,37 @@ +function bidsQA(opt, varargin) + % + % Run QA on a BIDS dataset. + % + + % (C) Copyright 2024 Remi Gau + + dsDesc = bids.util.jsondecode(fullfile(opt.dir.input, 'dataset_description.json')); + + switch dsDesc.DatasetType + + case 'raw' + + [BIDS, opt] = setUpWorkflow(opt, 'QA', opt.dir.input); + + bids.diagnostic(BIDS, ... + 'split_by', {'task'}, ... + 'output_path', fullfile(opt.dir.output, 'reports'), ... + 'verbose', false); + + case 'derivative' + + pipelines = {}; + for i = 1:numel(dsDesc.GeneratedBy) + pipelines {end+1} = dsDesc.GeneratedBy(i).Name; %#ok<*AGROW> + end + + if ismember('MRIQC', pipelines) + opt.dir.mriqc = opt.dir.input; + bidsQAmriqc(opt, 'T1w'); + bidsQAmriqc(opt, 'bold'); + + elseif any(ismember({'fMRIPrep', 'bidspm'}, pipelines)) + bidsQApreproc(opt, varargin{:}); + end + end +end diff --git a/src/workflows/bidsQApreproc.m b/src/workflows/QA/bidsQAbidspm.m similarity index 98% rename from src/workflows/bidsQApreproc.m rename to src/workflows/QA/bidsQAbidspm.m index 7d28144c5..e865d5f07 100644 --- a/src/workflows/bidsQApreproc.m +++ b/src/workflows/QA/bidsQAbidspm.m @@ -1,9 +1,10 @@ -function bidsQApreproc(opt) +function bidsQAbidspm(opt) % + % Run QA on preproc output generated by bidspm. % % USAGE:: % - % bidsQApreproc(opt) + % bidsQAbidspm(opt) % % :type opt: structure % :param opt: Options chosen for the analysis. @@ -13,6 +14,7 @@ function bidsQApreproc(opt) % For each run works on the realigned (and unwarped) data: % % - plots motion, global signal, framewise displacement + % % (C) Copyright 2020 bidspm developers diff --git a/src/QA/mriqcQA.m b/src/workflows/QA/bidsQAmriqc.m similarity index 100% rename from src/QA/mriqcQA.m rename to src/workflows/QA/bidsQAmriqc.m diff --git a/src/workflows/QA/bidsQApreproc.m b/src/workflows/QA/bidsQApreproc.m new file mode 100644 index 000000000..7212993b2 --- /dev/null +++ b/src/workflows/QA/bidsQApreproc.m @@ -0,0 +1,129 @@ +function bidsQApreproc(opt, varargin) + % + % Use confounds preprocessed datasets to estimate + % the number of timepoints with values superior to threshold. + % + % + % plots the proportion of timepoints per run + % (to identify runs with that goes above a limit) + + % (C) Copyright 2024 bidspm developers + + % TODO + % - handle opt.bidsFilterFile + % + % could also plots the sum of the proportion of timepoints over the 4 runs + % to identify participants whoe move a lot + % but for whom no run goes above the threshold + + args = inputParser; + addParameter(args, 'metric', 'framewise_displacement'); + addParameter(args, 'threshold', 0.2); + parse(args, varargin{:}); + metric = args.Results.metric; + threshold = args.Results.threshold; + + layoutFilter = struct('modality', {{'func'}}); + + [BIDS, opt] = setUpWorkflow(opt, 'QA', opt.dir.input, true, false, layoutFilter); + + if ~isfield(opt, 'taskName') + opt.taskName = '.*'; + end + + %% Load confounds to collect the metric for each run / subject + % stores only the metric of interest + + df = struct('sub', '', metric, []); + + for iSub = 1:numel(opt.subjects) + + subLabel = opt.subjects{iSub}; + + printProcessingSubject(iSub, subLabel, opt); + + filter = fileFilterForBold(opt, subLabel, 'confounds'); + + files = bids.query(BIDS, 'data', filter); + + for iRun = 1:numel(files) + + try + data = spm_load(files{iRun}); + catch + continue + end + + if ~ismember(metric, fieldnames(data)) + continue + end + + df(iSub).sub = subLabel; + df(iSub).(metric){iRun} = data.(metric); + + end + + end + + nbSubjects = numel(df); + + %% plot proportion datapoint with > threshold + for iSub = 1:nbSubjects + if isempty(df(iSub).(metric)) + continue + end + values = df(iSub).(metric); + df(iSub).proportion_censored = cellfun(@(x) sum(x > threshold) / sum(~isnan(x)), values); + end + + close all; + + plotFigure(df, metric, threshold); + +end + +function plotFigure(df, metric, threshold) + + nbSubjects = numel(df); + + figure('name', metric); + + hold on; + + x_tick_label = {df.sub}; + xtick = linspace(1.5, nbSubjects + 0.5, nbSubjects); + + % collect subjects average + avg = []; + for iSub = 1:nbSubjects + y = df(iSub).proportion_censored; + + avg(iSub) = mean(y); + if numel(y) == 1 + x = xtick(iSub); + elseif numel(y) == 2 + x = linspace(xtick(iSub) - 0.2, xtick(iSub) + 0.2, numel(y)); + elseif numel(y) == 3 + x = linspace(xtick(iSub) - 0.3, xtick(iSub) + 0.3, numel(y)); + else + x = linspace(xtick(iSub) - 0.35, xtick(iSub) + 0.35, numel(y)); + end + bar(x, y, 0.9); + plot(xtick(iSub), mean(y), 'k.', 'MarkerSize', 10); + end + + % plot group median + plot([1 nbSubjects + 1], [median(avg) median(avg)], '-r'); + + ylabel(sprintf('proportion time points %s > %0.5f / run', ... + strrep(metric, '_', ' '), ... + threshold)); + xlabel('subject'); + + set(gca, 'xtick', xtick, ... + 'xticklabel', x_tick_label, ... + 'ytick', 0:.05:1, ... + 'yticklabel', 0:.05:1, ... + 'fontsize', 8, 'tickdir', 'out'); + +end diff --git a/src/workflows/preproc/bidsResliceTpmToFunc.m b/src/workflows/preproc/bidsResliceTpmToFunc.m index fc16cfaa7..1e3531e18 100644 --- a/src/workflows/preproc/bidsResliceTpmToFunc.m +++ b/src/workflows/preproc/bidsResliceTpmToFunc.m @@ -16,7 +16,7 @@ function bidsResliceTpmToFunc(opt) % by :func:`bidsSpatialPrepro` % or :func:`bidsSegmentSkullStrip`. % - % It is necessary to run this workflow before running the :func:`bidsQApreproc` pipeline + % It is necessary to run this workflow before running the :func:`bidsQAbidspm` pipeline % as the computation of the tSNR by ``spmup`` requires the TPMs to have the same dimension % as the functional. % diff --git a/src/workflows/preproc/bidsSpatialPrepro.m b/src/workflows/preproc/bidsSpatialPrepro.m index 15be9d3d2..4be8d4253 100644 --- a/src/workflows/preproc/bidsSpatialPrepro.m +++ b/src/workflows/preproc/bidsSpatialPrepro.m @@ -155,7 +155,7 @@ transferMetadataFromJson(createdFiles); - bidsQApreproc(opt); + bidsQAbidspm(opt); end diff --git a/tests/tests_workflows/QA/test_bidsQA.m b/tests/tests_workflows/QA/test_bidsQA.m new file mode 100644 index 000000000..e53d749ed --- /dev/null +++ b/tests/tests_workflows/QA/test_bidsQA.m @@ -0,0 +1,80 @@ +function test_suite = test_bidsQA %#ok<*STOUT> + + % (C) Copyright 2024 bidspm developers + + try % assignment of 'localfunctions' is necessary in Matlab >= 2016 + test_functions = localfunctions(); %#ok<*NASGU> + catch % no problem; early Matlab versions can use initTestSuite fine + end + initTestSuite; +end + +function test_bidsQA_raw() + + % opt = setOptions('vislocalizer'); + % + % expectedOutput = fullfile(opt.dir.output, ... + % 'reports', ... + % 'bidspm-raw_split_by-task.png'); + % + % if exist(expectedOutput, 'file') + % delete(expectedOutput) + % end + % + % bidsQA(opt) + % + % assertEqual(exist(expectedOutput, 'file'), 2) + +end + +function test_bidsQA_mriqc() + + ds000114mriqc = spm_file(fullfile(getTestDir(), ... + '..', 'demos', 'openneuro', ... + 'inputs', 'ds000114-mriqc'), 'cpath'); + + opt.dir.input = ds000114mriqc; + + bidsQA(opt); + +end + +function test_bidsQA_bidspm() + + % opt = setOptions('vislocalizer'); + % + % opt.dir.input = opt.dir.preproc; + % + % opt.dir.output = ''; + % opt.dir.derivatives = tempName(); + % + % opt = checkOptions(opt); + % + % threshold = 0.0001; + % metric = 'rot_x'; + % + % bidsQA(opt, 'metric', metric, 'threshold', threshold); + +end + +function test_bidsQA_fmriprep() + + ds000114fmriprep = spm_file(fullfile(getTestDir(), ... + '..', 'demos', 'openneuro', ... + 'inputs', 'ds000114-fmriprep'), 'cpath'); + + opt.dir.input = ds000114fmriprep; + opt.dir.preproc = ds000114fmriprep; + opt.dir.derivatives = tempName(); + opt.pipeline.type = 'preproc'; + opt.verbosity = 0; + + opt.taskName = {'overtverbgeneration', 'linebisection'}; + + opt = checkOptions(opt); + + opt.subjects = {'01', '02', '03', '10'}; + + bidsQA(opt); + +end diff --git a/tests/tests_QA/test_mriqcQA.m b/tests/tests_workflows/QA/test_bidsQAmriqc.m similarity index 81% rename from tests/tests_QA/test_mriqcQA.m rename to tests/tests_workflows/QA/test_bidsQAmriqc.m index b5daa4821..9846a736b 100644 --- a/tests/tests_QA/test_mriqcQA.m +++ b/tests/tests_workflows/QA/test_bidsQAmriqc.m @@ -1,4 +1,4 @@ -function test_suite = test_mriqcQA %#ok<*STOUT> +function test_suite = test_bidsQAmriqc %#ok<*STOUT> % (C) Copyright 2022 bidspm developers try % assignment of 'localfunctions' is necessary in Matlab >= 2016 test_functions = localfunctions(); %#ok<*NASGU> @@ -11,6 +11,6 @@ function test_mriqc_basic() opt.dir.mriqc = fullfile(getTestDataDir, 'tsv_files'); - mriqcQA(opt, 'T1w'); + bidsQAmriqc(opt, 'T1w'); end From 87ab6bc31ff0ba49f1cf2af6bdf74de62b150f99 Mon Sep 17 00:00:00 2001 From: Remi Gau Date: Wed, 31 Jul 2024 12:17:26 +0200 Subject: [PATCH 2/7] start adding QA workflow --- src/workflows/QA/bidsQApreproc.m | 1 + 1 file changed, 1 insertion(+) diff --git a/src/workflows/QA/bidsQApreproc.m b/src/workflows/QA/bidsQApreproc.m index 7212993b2..c89ab61bf 100644 --- a/src/workflows/QA/bidsQApreproc.m +++ b/src/workflows/QA/bidsQApreproc.m @@ -11,6 +11,7 @@ function bidsQApreproc(opt, varargin) % TODO % - handle opt.bidsFilterFile + % - save collected values to tsv to make it easier to apply new threshold % % could also plots the sum of the proportion of timepoints over the 4 runs % to identify participants whoe move a lot From 883cb05a6b2dc70a452b6415267b95d7e20e42ec Mon Sep 17 00:00:00 2001 From: Remi Gau Date: Wed, 31 Jul 2024 12:20:21 +0200 Subject: [PATCH 3/7] update changelog --- CHANGELOG.md | 1 + 1 file changed, 1 insertion(+) diff --git a/CHANGELOG.md b/CHANGELOG.md index 3e526c589..55087afe8 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -24,6 +24,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added +* [ENH] add QA workflow `bidsQA` to find outliers in MRIQC output or to view number of outlier timepoints (for a given metric and threshold) in each functional run #1297 by @Remi-Gau * [ENH] add support for one-way ANOVA across groups at the group level #1296 by @Remi-Gau * [ENH] allow for 2 sample T-Test, within group T-Test and one-way ANOVA to ne more flexible with respect to what praticipants.tsv column to use to allocate subjects in each group #1296 by @Remi-Gau * [ENH] make `addConfoundsToDesignMatrix` a method of `BidsModel` #1294 by @Remi-Gau From 65c70c2118c15bd3c480531c2218ba6d79308f1c Mon Sep 17 00:00:00 2001 From: Remi Gau Date: Wed, 31 Jul 2024 12:44:12 +0200 Subject: [PATCH 4/7] move files to slow --- .github/workflows/tests.yml | 1 + .github/workflows/tests_octave.yml | 6 ++ docs/source/API.rst | 13 ++- src/workflows/QA/bidsQA.m | 19 +++- src/workflows/QA/bidsQAmriqc.m | 8 +- src/workflows/QA/bidsQApreproc.m | 36 ++++++- .../tests_workflows/QA/test_bidsQA.m | 100 ++++++++++++++++++ tests/tests_workflows/QA/test_bidsQA.m | 80 -------------- tests/tests_workflows/QA/test_bidsQAmriqc.m | 16 --- 9 files changed, 171 insertions(+), 108 deletions(-) create mode 100644 tests/tests_slow/tests_workflows/QA/test_bidsQA.m delete mode 100644 tests/tests_workflows/QA/test_bidsQA.m delete mode 100644 tests/tests_workflows/QA/test_bidsQAmriqc.m diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index 818bac6c6..804e021ab 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -121,6 +121,7 @@ jobs: run: | cd demos/openneuro/ make data_ds000114_mriqc + make data_ds000114_fmriprep - name: Prepare test data unix run: | diff --git a/.github/workflows/tests_octave.yml b/.github/workflows/tests_octave.yml index 92a46c57f..0de29b257 100644 --- a/.github/workflows/tests_octave.yml +++ b/.github/workflows/tests_octave.yml @@ -94,6 +94,12 @@ jobs: unzip download mv moae_fmriprep fmriprep + - name: Get data for testing QA + run: | + cd demos/openneuro/ + make data_ds000114_mriqc + make data_ds000114_fmriprep + - name: Prepare test data run: | cd tests diff --git a/docs/source/API.rst b/docs/source/API.rst index 4f16e8e37..d226f4f4e 100644 --- a/docs/source/API.rst +++ b/docs/source/API.rst @@ -33,7 +33,6 @@ workflows .. autofunction:: src.workflows.bidsCheckVoxelSize .. autofunction:: src.workflows.bidsCopyInputFolder .. autofunction:: src.workflows.bidsInverseNormalize -.. autofunction:: src.workflows.bidsQAbidspm .. autofunction:: src.workflows.bidsRename .. autofunction:: src.workflows.bidsReport @@ -42,6 +41,13 @@ roi .. autofunction:: src.workflows.roi.bidsCreateROI .. autofunction:: src.workflows.roi.bidsRoiBasedGLM +QA +-- +.. autofunction:: src.workflows.QA.bidsQA +.. autofunction:: src.workflows.QA.bidsQAbidspm +.. autofunction:: src.workflows.QA.bidsQAmriqc +.. autofunction:: src.workflows.QA.bidsQApreproc + lesion ------ .. autofunction:: src.workflows.lesion.bidsLesionAbnormalitiesDetection @@ -103,7 +109,6 @@ stats .. autofunction:: src.batches.stats.setBatchSubjectLevelContrasts .. autofunction:: src.batches.stats.setBatchSubjectLevelGLMSpec .. autofunction:: src.batches.stats.setBatchSubjectLevelResults -.. autofunction:: src.batches.stats.setBatchTwoSampleTTest preproc ------- @@ -157,7 +162,6 @@ QA .. autofunction:: src.QA.computeFDandRMS .. autofunction:: src.QA.computeRobustOutliers .. autofunction:: src.QA.createDesignMatrix -.. autofunction:: src.QA.mriqcQA .. autofunction:: src.QA.plotConfounds .. autofunction:: src.QA.plotEvents .. autofunction:: src.QA.plotRoiTimeCourse @@ -204,6 +208,9 @@ bidspm data ---- +__pycache__ +----------- + cli === .. autofunction:: src.cli.baseInputParser diff --git a/src/workflows/QA/bidsQA.m b/src/workflows/QA/bidsQA.m index 13e283b35..79215282f 100644 --- a/src/workflows/QA/bidsQA.m +++ b/src/workflows/QA/bidsQA.m @@ -1,10 +1,25 @@ -function bidsQA(opt, varargin) +function filename = bidsQA(opt, varargin) % % Run QA on a BIDS dataset. % + % Find outliers in MRIQC output + % or to view number of outlier timepoints + % (for a given metric and threshold) in each functional run + % + % + % USAGE:: + % + % figurePath = bidsQA(opt, varargin); + % + % :param opt: Options chosen for the analysis. + % See :func:`checkOptions`. + % :type opt: structure + % % (C) Copyright 2024 Remi Gau + filename = ''; + dsDesc = bids.util.jsondecode(fullfile(opt.dir.input, 'dataset_description.json')); switch dsDesc.DatasetType @@ -31,7 +46,7 @@ function bidsQA(opt, varargin) bidsQAmriqc(opt, 'bold'); elseif any(ismember({'fMRIPrep', 'bidspm'}, pipelines)) - bidsQApreproc(opt, varargin{:}); + filename = bidsQApreproc(opt, varargin{:}); end end end diff --git a/src/workflows/QA/bidsQAmriqc.m b/src/workflows/QA/bidsQAmriqc.m index 738b3b3d2..dc3c3084b 100644 --- a/src/workflows/QA/bidsQAmriqc.m +++ b/src/workflows/QA/bidsQAmriqc.m @@ -1,4 +1,4 @@ -function mriqcQA(opt, suffix) +function bidsQAmriqc(opt, suffix) % % uses the report from MRIQC (bold and T1) identify outliers using % robust statistics (interquartile range). @@ -8,7 +8,7 @@ function mriqcQA(opt, suffix) % % USAGE:: % - % mriqcQA(opt, suffix); + % bidsQAmriqc(opt, suffix); % % :type opt: structure % :param opt: Options chosen for the analysis. @@ -18,8 +18,8 @@ function mriqcQA(opt, suffix) % % opt.dir.mriqc = '/home/remi/gin/dataset/derivatives/mriqc'; % - % mriqcQA(opt, 'T1w'); - % mriqcQA(opt, 'bold'); + % bidsQAmriqc(opt, 'T1w'); + % bidsQAmriqc(opt, 'bold'); % % Dependencies (in case you want to use it as standalone): % diff --git a/src/workflows/QA/bidsQApreproc.m b/src/workflows/QA/bidsQApreproc.m index c89ab61bf..d3b012c3c 100644 --- a/src/workflows/QA/bidsQApreproc.m +++ b/src/workflows/QA/bidsQApreproc.m @@ -1,11 +1,25 @@ -function bidsQApreproc(opt, varargin) +function filename = bidsQApreproc(opt, varargin) % % Use confounds preprocessed datasets to estimate % the number of timepoints with values superior to threshold. % + % Plots the proportion of timepoints per run + % (to identify runs with that goes above a limit). + % + % USAGE:: + % + % figurePath = bidsQApreproc(opt, 'metric', metric, 'threshold', threshold); + % + % :param opt: Options chosen for the analysis. + % See :func:`checkOptions`. + % :type opt: structure + % + % :param metric: default = 'framewise_displacement' + % :type metric: char + % + % :param threshold: default = 0.2 + % :param threshold: numeric % - % plots the proportion of timepoints per run - % (to identify runs with that goes above a limit) % (C) Copyright 2024 bidspm developers @@ -81,6 +95,9 @@ function bidsQApreproc(opt, varargin) plotFigure(df, metric, threshold); + output_path = fullfile(opt.dir.output, 'reports'); + filename = printFigure(output_path, metric); + end function plotFigure(df, metric, threshold) @@ -128,3 +145,16 @@ function plotFigure(df, metric, threshold) 'fontsize', 8, 'tickdir', 'out'); end + +function filename = printFigure(outputPath, metric) + bids.util.mkdir(outputPath); + + filename = regexprep(['outliers_', metric, '.png'], ' - ', '_'); + filename = regexprep(filename, ' ', '-'); + filename = regexprep(filename, '[\(\)]', ''); + filename = fullfile(outputPath, filename); + + print(filename, '-dpng'); + fprintf('Figure saved:\n\t%s\n', filename); + +end diff --git a/tests/tests_slow/tests_workflows/QA/test_bidsQA.m b/tests/tests_slow/tests_workflows/QA/test_bidsQA.m new file mode 100644 index 000000000..7dabedc20 --- /dev/null +++ b/tests/tests_slow/tests_workflows/QA/test_bidsQA.m @@ -0,0 +1,100 @@ +function test_suite = test_bidsQA %#ok<*STOUT> + + % (C) Copyright 2024 bidspm developers + + try % assignment of 'localfunctions' is necessary in Matlab >= 2016 + test_functions = localfunctions(); %#ok<*NASGU> + catch % no problem; early Matlab versions can use initTestSuite fine + end + initTestSuite; +end + +function test_bidsQA_raw() + + % markTestAs('slow'); + % + % opt = setOptions('vislocalizer'); + % + % expectedOutput = fullfile(opt.dir.output, ... + % 'reports', ... + % 'bidspm-raw_split_by-task.png'); + % + % if exist(expectedOutput, 'file') + % delete(expectedOutput) + % end + % + % bidsQA(opt) + % + % assertEqual(exist(expectedOutput, 'file'), 2) + +end + +function test_bidsQA_mriqc() + + % markTestAs('slow'); + % + % if ispc() + % moxunit_throw_test_skipped_exception('requires datalad setup'); + % end + % + % ds000114mriqc = spm_file(fullfile(getTestDir(), ... + % '..', 'demos', 'openneuro', ... + % 'inputs', 'ds000114-mriqc'), 'cpath'); + % + % opt.dir.input = ds000114mriqc; + % + % bidsQA(opt); + +end + +function test_bidsQA_bidspm() + + markTestAs('slow'); + + if ispc() + moxunit_throw_test_skipped_exception('requires datalad setup'); + end + + opt = setOptions('vislocalizer'); + + opt.dir.input = opt.dir.preproc; + + opt.dir.output = ''; + opt.dir.derivatives = tempName(); + + opt = checkOptions(opt); + + threshold = 0.0001; + metric = 'rot_x'; + + filename = bidsQA(opt, 'metric', metric, 'threshold', threshold); + + assertEqual(exist(filename, 'file'), 2); + +end + +function test_bidsQA_fmriprep() + + markTestAs('slow'); + + ds000114fmriprep = spm_file(fullfile(getTestDir(), ... + '..', 'demos', 'openneuro', ... + 'inputs', 'ds000114-fmriprep'), 'cpath'); + + opt.dir.input = ds000114fmriprep; + opt.dir.preproc = ds000114fmriprep; + opt.dir.derivatives = tempName(); + opt.pipeline.type = 'preproc'; + opt.verbosity = 0; + + opt.taskName = {'overtverbgeneration', 'linebisection'}; + + opt = checkOptions(opt); + + opt.subjects = {'01', '02', '03', '10'}; + + filename = bidsQA(opt); + + assertEqual(exist(filename, 'file'), 2); + +end diff --git a/tests/tests_workflows/QA/test_bidsQA.m b/tests/tests_workflows/QA/test_bidsQA.m deleted file mode 100644 index e53d749ed..000000000 --- a/tests/tests_workflows/QA/test_bidsQA.m +++ /dev/null @@ -1,80 +0,0 @@ -function test_suite = test_bidsQA %#ok<*STOUT> - - % (C) Copyright 2024 bidspm developers - - try % assignment of 'localfunctions' is necessary in Matlab >= 2016 - test_functions = localfunctions(); %#ok<*NASGU> - catch % no problem; early Matlab versions can use initTestSuite fine - end - initTestSuite; -end - -function test_bidsQA_raw() - - % opt = setOptions('vislocalizer'); - % - % expectedOutput = fullfile(opt.dir.output, ... - % 'reports', ... - % 'bidspm-raw_split_by-task.png'); - % - % if exist(expectedOutput, 'file') - % delete(expectedOutput) - % end - % - % bidsQA(opt) - % - % assertEqual(exist(expectedOutput, 'file'), 2) - -end - -function test_bidsQA_mriqc() - - ds000114mriqc = spm_file(fullfile(getTestDir(), ... - '..', 'demos', 'openneuro', ... - 'inputs', 'ds000114-mriqc'), 'cpath'); - - opt.dir.input = ds000114mriqc; - - bidsQA(opt); - -end - -function test_bidsQA_bidspm() - - % opt = setOptions('vislocalizer'); - % - % opt.dir.input = opt.dir.preproc; - % - % opt.dir.output = ''; - % opt.dir.derivatives = tempName(); - % - % opt = checkOptions(opt); - % - % threshold = 0.0001; - % metric = 'rot_x'; - % - % bidsQA(opt, 'metric', metric, 'threshold', threshold); - -end - -function test_bidsQA_fmriprep() - - ds000114fmriprep = spm_file(fullfile(getTestDir(), ... - '..', 'demos', 'openneuro', ... - 'inputs', 'ds000114-fmriprep'), 'cpath'); - - opt.dir.input = ds000114fmriprep; - opt.dir.preproc = ds000114fmriprep; - opt.dir.derivatives = tempName(); - opt.pipeline.type = 'preproc'; - opt.verbosity = 0; - - opt.taskName = {'overtverbgeneration', 'linebisection'}; - - opt = checkOptions(opt); - - opt.subjects = {'01', '02', '03', '10'}; - - bidsQA(opt); - -end diff --git a/tests/tests_workflows/QA/test_bidsQAmriqc.m b/tests/tests_workflows/QA/test_bidsQAmriqc.m deleted file mode 100644 index 9846a736b..000000000 --- a/tests/tests_workflows/QA/test_bidsQAmriqc.m +++ /dev/null @@ -1,16 +0,0 @@ -function test_suite = test_bidsQAmriqc %#ok<*STOUT> - % (C) Copyright 2022 bidspm developers - try % assignment of 'localfunctions' is necessary in Matlab >= 2016 - test_functions = localfunctions(); %#ok<*NASGU> - catch % no problem; early Matlab versions can use initTestSuite fine - end - initTestSuite; -end - -function test_mriqc_basic() - - opt.dir.mriqc = fullfile(getTestDataDir, 'tsv_files'); - - bidsQAmriqc(opt, 'T1w'); - -end From c52624a279bfa1f8e19f0d502683feca72a79338 Mon Sep 17 00:00:00 2001 From: Remi Gau Date: Wed, 31 Jul 2024 12:55:48 +0200 Subject: [PATCH 5/7] activate codecov tokens --- .github/workflows/tests.yml | 2 +- .github/workflows/tests_octave.yml | 2 +- .github/workflows/tests_windows.yml | 2 +- demos/openneuro/ds000114_preproc_run.m | 12 ++++++------ 4 files changed, 9 insertions(+), 9 deletions(-) diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index 804e021ab..e36d265ee 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -162,7 +162,7 @@ jobs: flags: ${{ matrix.os }}_matlab-${{ matrix.matlab }}_${{ matrix.mode }} name: codecov-matlab fail_ci_if_error: false - # token: ${{ secrets.CODECOV_TOKEN }} # not required but might help API rate limits + token: ${{ secrets.CODECOV_TOKEN }} - name: Run system tests MATLAB ${{ matrix.script }} if: matrix.test_type == 'system' diff --git a/.github/workflows/tests_octave.yml b/.github/workflows/tests_octave.yml index 0de29b257..db53598b2 100644 --- a/.github/workflows/tests_octave.yml +++ b/.github/workflows/tests_octave.yml @@ -155,7 +155,7 @@ jobs: flags: octave name: codecov-octave_${{ matrix.mode }} fail_ci_if_error: false - # token: ${{ secrets.CODECOV_TOKEN }} # not required but might help API rate limits + # token: ${{ secrets.CODECOV_TOKEN }} - name: Run system tests octave ${{ matrix.script }} if: matrix.test_type == 'system' diff --git a/.github/workflows/tests_windows.yml b/.github/workflows/tests_windows.yml index f14f2ebce..10a57409c 100644 --- a/.github/workflows/tests_windows.yml +++ b/.github/workflows/tests_windows.yml @@ -131,7 +131,7 @@ jobs: flags: ${{ matrix.os }}_matlab-${{ matrix.matlab }}_${{ matrix.mode }} name: codecov-matlab fail_ci_if_error: false - # token: ${{ secrets.CODECOV_TOKEN }} # not required but might help API rate limits + # token: ${{ secrets.CODECOV_TOKEN }} - name: Run system tests MATLAB ${{ matrix.script }} if: matrix.test_type == 'system' diff --git a/demos/openneuro/ds000114_preproc_run.m b/demos/openneuro/ds000114_preproc_run.m index 7e6be75c4..e9d299069 100644 --- a/demos/openneuro/ds000114_preproc_run.m +++ b/demos/openneuro/ds000114_preproc_run.m @@ -10,7 +10,7 @@ participant_label = {'04', '05'}; TASK = 'linebisection'; -session_to_select = {'test', 'retest'}; +session_to_select = {'test'}; root_dir = fileparts(mfilename('fullpath')); bids_dir = fullfile(root_dir, 'inputs', 'ds000114'); @@ -19,11 +19,11 @@ %% Copy % To make sure we still have copied all the data -% bidspm(bids_dir, output_dir, 'subject', ... -% 'participant_label', participant_label, ... -% 'action', 'copy', ... -% 'skip_validation', true, ... -% 'verbosity', 3); +bidspm(bids_dir, output_dir, 'subject', ... + 'participant_label', participant_label, ... + 'action', 'copy', ... + 'skip_validation', true, ... + 'verbosity', 3); %% Set filter and preprocess bids_filter_file = struct('bold', struct('modality', 'func', ... From 3d6544a0da9fb713e5cb38ff7eb677b7d08d0040 Mon Sep 17 00:00:00 2001 From: Remi Gau Date: Wed, 31 Jul 2024 12:58:32 +0200 Subject: [PATCH 6/7] update doc --- README.md | 2 +- docs/source/QA.rst | 7 ++++--- 2 files changed, 5 insertions(+), 4 deletions(-) diff --git a/README.md b/README.md index 4019766d6..92985a4fb 100644 --- a/README.md +++ b/README.md @@ -153,7 +153,7 @@ then bidspm has some automated workflows to perform amongst other things: All (well almost all) preprocessed outputs are saved as BIDS derivatives with BIDS compliant filenames. -### Quality control: +### Quality control - anatomical data (work in progress) - functional data (work in progress) diff --git a/docs/source/QA.rst b/docs/source/QA.rst index 6cb0b016c..cb649b0f9 100644 --- a/docs/source/QA.rst +++ b/docs/source/QA.rst @@ -5,11 +5,11 @@ Quality control The illustrations in this section mix what the files created by each workflow and the functions and are called by it. - In this sense they are not pure DAGs (directed acyclic graphs) as the ``*.m`` files - mentioned in them already exist. + In this sense they are not pure DAGs (directed acyclic graphs) + as the ``*.m`` files mentioned in them already exist. -- :func:`anatQA` +- :func:`bidsQA` - :func:`bidsQAbidspm` .. _fig_spatialPrepro-reports: @@ -18,5 +18,6 @@ Quality control workflows for QA as part of the spatial preprocessing workflow +- :func:`anatQA` - :func:`computeDesignEfficiency` - :func:`plotEvents` From 9933b6c0a19e0f32680cbde518e62478dba9a759 Mon Sep 17 00:00:00 2001 From: Remi Gau Date: Wed, 31 Jul 2024 13:21:41 +0200 Subject: [PATCH 7/7] fix set up --- .github/workflows/tests.yml | 1 + .github/workflows/tests_octave.yml | 7 +- .../tests_workflows/QA/test_bidsQA.m | 64 +++++++++---------- 3 files changed, 39 insertions(+), 33 deletions(-) diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index e36d265ee..a69349405 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -118,6 +118,7 @@ jobs: mv moae_fmriprep fmriprep - name: Get data for testing QA + if: matrix.os == 'ubuntu-latest' run: | cd demos/openneuro/ make data_ds000114_mriqc diff --git a/.github/workflows/tests_octave.yml b/.github/workflows/tests_octave.yml index db53598b2..93bd47e2b 100644 --- a/.github/workflows/tests_octave.yml +++ b/.github/workflows/tests_octave.yml @@ -61,7 +61,7 @@ jobs: - name: Install dependencies run: | sudo apt-get -y -qq update - sudo apt-get -y install unzip wget + sudo apt-get -y install unzip wget git-annex - name: Install Node uses: actions/setup-node@v4 @@ -79,6 +79,11 @@ jobs: submodules: recursive fetch-depth: 0 + - name: Install datalad + run: | + python -m pip install --upgrade pip setuptools + pip install datalad + - name: Install validators run: make install diff --git a/tests/tests_slow/tests_workflows/QA/test_bidsQA.m b/tests/tests_slow/tests_workflows/QA/test_bidsQA.m index 7dabedc20..829a2494c 100644 --- a/tests/tests_slow/tests_workflows/QA/test_bidsQA.m +++ b/tests/tests_slow/tests_workflows/QA/test_bidsQA.m @@ -11,50 +11,46 @@ function test_bidsQA_raw() - % markTestAs('slow'); - % - % opt = setOptions('vislocalizer'); - % - % expectedOutput = fullfile(opt.dir.output, ... - % 'reports', ... - % 'bidspm-raw_split_by-task.png'); - % - % if exist(expectedOutput, 'file') - % delete(expectedOutput) - % end - % - % bidsQA(opt) - % - % assertEqual(exist(expectedOutput, 'file'), 2) + markTestAs('slow'); -end + opt = setOptions('vislocalizer'); -function test_bidsQA_mriqc() + expectedOutput = fullfile(opt.dir.output, ... + 'reports', ... + 'bidspm-raw_split_by-task.png'); + + if exist(expectedOutput, 'file') + delete(expectedOutput); + end - % markTestAs('slow'); - % - % if ispc() - % moxunit_throw_test_skipped_exception('requires datalad setup'); - % end - % - % ds000114mriqc = spm_file(fullfile(getTestDir(), ... - % '..', 'demos', 'openneuro', ... - % 'inputs', 'ds000114-mriqc'), 'cpath'); - % - % opt.dir.input = ds000114mriqc; - % - % bidsQA(opt); + bidsQA(opt); + + assertEqual(exist(expectedOutput, 'file'), 2); end -function test_bidsQA_bidspm() +function test_bidsQA_mriqc() markTestAs('slow'); - if ispc() + if ispc() || ismac() moxunit_throw_test_skipped_exception('requires datalad setup'); end + ds000114mriqc = spm_file(fullfile(getTestDir(), ... + '..', 'demos', 'openneuro', ... + 'inputs', 'ds000114-mriqc'), 'cpath'); + + opt.dir.input = ds000114mriqc; + + bidsQA(opt); + +end + +function test_bidsQA_bidspm() + + markTestAs('slow'); + opt = setOptions('vislocalizer'); opt.dir.input = opt.dir.preproc; @@ -77,6 +73,10 @@ function test_bidsQA_fmriprep() markTestAs('slow'); + if ispc() || ismac() + moxunit_throw_test_skipped_exception('requires datalad setup'); + end + ds000114fmriprep = spm_file(fullfile(getTestDir(), ... '..', 'demos', 'openneuro', ... 'inputs', 'ds000114-fmriprep'), 'cpath');