diff --git a/Dockerfile b/Dockerfile index f236bea..bcd1e2f 100644 --- a/Dockerfile +++ b/Dockerfile @@ -237,7 +237,7 @@ ENV IS_DOCKER_8395080871=1 RUN ldconfig WORKDIR /tmp -ENTRYPOINT ["/opt/conda/envs/fmripost_aroma/bin/fmripost_aroma"] +ENTRYPOINT ["/opt/conda/envs/fmripost_aroma/bin/fmripost-aroma"] ARG BUILD_DATE ARG VCS_REF diff --git a/docs/spaces.rst b/docs/spaces.rst index 1c90ea3..cb87146 100644 --- a/docs/spaces.rst +++ b/docs/spaces.rst @@ -37,7 +37,7 @@ stereotactic coordinate system. The most extended standard space for fMRI analyses is generally referred to MNI. For instance, to instruct *fMRIPost-AROMA* to use the MNI template brain distributed with FSL as coordinate reference the option will read as follows: ``--output-spaces MNI152NLin6Asym``. -By default, *fMRIPost-AROMA* uses ``MNI152NLin2009cAsym`` as spatial-standardization reference. +By default, *fMRIPost-AROMA* uses ``MNI152NLin6Asym`` as spatial-standardization reference. Valid template identifiers (``MNI152NLin6Asym``, ``MNI152NLin2009cAsym``, etc.) come from the `TemplateFlow repository `__. diff --git a/pyproject.toml b/pyproject.toml index fae85c8..5c147e4 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -59,7 +59,7 @@ Issues = "https://github.com/nipreps/fmripost-aroma/issues" Source = "https://github.com/nipreps/fmripost-aroma" [project.scripts] -fmripost-aroma = "fmripost_aroma.cli:aroma" +fmripost-aroma = "fmripost_aroma.cli.run:main" [tool.hatch.metadata] allow-direct-references = true diff --git a/src/fmripost_aroma/cli/parser.py b/src/fmripost_aroma/cli/parser.py index 3542ee8..92c5285 100644 --- a/src/fmripost_aroma/cli/parser.py +++ b/src/fmripost_aroma/cli/parser.py @@ -195,7 +195,6 @@ def _bids_filter(value, parser): '--derivatives', action=ToDict, metavar='PACKAGE=PATH', - type=str, nargs='+', help=( 'Search PATH(s) for pre-computed derivatives. ' @@ -332,6 +331,17 @@ def _bids_filter(value, parser): default=False, help='Skip generation of HTML and LaTeX formatted citation with pandoc', ) + g_outputs.add_argument( + '--aggregate-session-reports', + dest='aggr_ses_reports', + action='store', + type=PositiveInt, + default=4, + help=( + "Maximum number of sessions aggregated in one subject's visual report. " + 'If exceeded, visual reports are split by session.' + ), + ) g_aroma = parser.add_argument_group('Options for running ICA_AROMA') g_aroma.add_argument( @@ -347,7 +357,7 @@ def _bids_filter(value, parser): ) g_aroma.add_argument( '--error-on-warnings', - dest='error_on_aroma_warnings', + dest='err_on_warn', action='store_true', default=False, help=( @@ -361,6 +371,7 @@ def _bids_filter(value, parser): nargs='+', choices=['aggr', 'nonaggr', 'orthaggr'], default=None, + dest='denoise_method', help='Denoising method to apply, if any.', ) @@ -530,6 +541,9 @@ def parse_args(args=None, namespace=None): work_dir = config.execution.work_dir version = config.environment.version + if config.execution.fmripost_aroma_dir is None: + config.execution.fmripost_aroma_dir = output_dir + # Wipe out existing work_dir if opts.clean_workdir and work_dir.exists(): from niworkflows.utils.misc import clean_directory @@ -569,7 +583,7 @@ def parse_args(args=None, namespace=None): validate_input_dir(config.environment.exec_env, opts.bids_dir, opts.participant_label) # Setup directories - config.execution.log_dir = config.execution.fmriprep_dir / 'logs' + config.execution.log_dir = config.execution.fmripost_aroma_dir / 'logs' # Check and create output and working directories config.execution.log_dir.mkdir(exist_ok=True, parents=True) work_dir.mkdir(exist_ok=True, parents=True) diff --git a/src/fmripost_aroma/cli/run.py b/src/fmripost_aroma/cli/run.py index 0dac556..eabb18a 100644 --- a/src/fmripost_aroma/cli/run.py +++ b/src/fmripost_aroma/cli/run.py @@ -204,8 +204,6 @@ def main(): ) errno = 0 finally: - from pkg_resources import resource_filename as pkgrf - # Code Carbon if config.execution.track_carbon: emissions: float = tracker.stop() @@ -217,11 +215,9 @@ def main(): # Generate reports phase failed_reports = generate_reports( - config.execution.participant_label, - config.execution.fmripost_aroma_dir, - config.execution.run_uuid, - config=pkgrf('fmripost_aroma', 'data/reports-spec.yml'), - packagename='fmripost_aroma', + subject_list=config.execution.participant_label, + output_dir=config.execution.fmripost_aroma_dir, + run_uuid=config.execution.run_uuid, ) write_derivative_description( config.execution.bids_dir, config.execution.fmripost_aroma_dir @@ -233,7 +229,7 @@ def main(): 'Report generation failed for %d subjects' % failed_reports, level='error', ) - sys.exit(int((errno + failed_reports) > 0)) + sys.exit(int((errno + len(failed_reports)) > 0)) def migas_exit() -> None: diff --git a/src/fmripost_aroma/cli/workflow.py b/src/fmripost_aroma/cli/workflow.py index c1c7bf7..9e153f2 100644 --- a/src/fmripost_aroma/cli/workflow.py +++ b/src/fmripost_aroma/cli/workflow.py @@ -36,14 +36,13 @@ def build_workflow(config_file, retval): """Create the Nipype Workflow that supports the whole execution graph.""" from pathlib import Path + from fmriprep.utils.bids import check_pipeline_version + from fmriprep.utils.misc import check_deps from niworkflows.utils.bids import collect_participants - from niworkflows.utils.misc import check_valid_fs_license from pkg_resources import resource_filename as pkgrf from fmripost_aroma import config from fmripost_aroma.reports.core import generate_reports - from fmripost_aroma.utils.bids import check_pipeline_version - from fmripost_aroma.utils.misc import check_deps from fmripost_aroma.workflows.base import init_fmripost_aroma_wf config.load(config_file) @@ -92,11 +91,9 @@ def build_workflow(config_file, retval): if config.execution.reports_only: build_log.log(25, 'Running --reports-only on participants %s', ', '.join(subject_list)) retval['return_code'] = generate_reports( - config.execution.participant_label, - config.execution.fmripost_aroma_dir, - config.execution.run_uuid, - config=pkgrf('fmripost_aroma', 'data/reports-spec.yml'), - packagename='fmripost_aroma', + subject_list=config.execution.participant_label, + output_dir=config.execution.fmripost_aroma_dir, + run_uuid=config.execution.run_uuid, ) return retval @@ -112,37 +109,10 @@ def build_workflow(config_file, retval): if config.execution.derivatives: init_msg += [f'Searching for derivatives: {config.execution.derivatives}.'] - if config.execution.fs_subjects_dir: - init_msg += [f"Pre-run FreeSurfer's SUBJECTS_DIR: {config.execution.fs_subjects_dir}."] - build_log.log(25, f"\n{' ' * 11}* ".join(init_msg)) retval['workflow'] = init_fmripost_aroma_wf() - # Check for FS license after building the workflow - if not check_valid_fs_license(): - from fmripost_aroma.utils.misc import fips_enabled - - if fips_enabled(): - build_log.critical( - """\ -ERROR: Federal Information Processing Standard (FIPS) mode is enabled on your system. \ -FreeSurfer (and thus fMRIPost-AROMA) cannot be used in FIPS mode. \ -Contact your system administrator for assistance.""" - ) - else: - build_log.critical( - """\ -ERROR: a valid license file is required for FreeSurfer to run. \ -fMRIPost-AROMA looked for an existing license file at several paths, in this order: \ -1) command line argument ``--fs-license-file``; \ -2) ``$FS_LICENSE`` environment variable; and \ -3) the ``$FREESURFER_HOME/license.txt`` path. \ -Get it (for free) by registering at https://surfer.nmr.mgh.harvard.edu/registration.html""" - ) - retval['return_code'] = 126 # 126 == Command invoked cannot execute. - return retval - # Check workflow for missing commands missing = check_deps(retval['workflow']) if missing: diff --git a/src/fmripost_aroma/config.py b/src/fmripost_aroma/config.py index b02d3b6..0aa49f1 100644 --- a/src/fmripost_aroma/config.py +++ b/src/fmripost_aroma/config.py @@ -225,6 +225,7 @@ def load(cls, settings, init=True, ignore=None): for k, v in settings.items(): if k in ignore or v is None: continue + if k in cls._paths: if isinstance(v, list | tuple): setattr(cls, k, [Path(val).absolute() for val in v]) @@ -232,6 +233,7 @@ def load(cls, settings, init=True, ignore=None): setattr(cls, k, {key: Path(val).absolute() for key, val in v.items()}) else: setattr(cls, k, Path(v).absolute()) + elif hasattr(cls, k): setattr(cls, k, v) @@ -250,17 +252,24 @@ def get(cls): for k, v in cls.__dict__.items(): if k.startswith('_') or v is None: continue + if callable(getattr(cls, k)): continue + if k in cls._paths: if isinstance(v, list | tuple): v = [str(val) for val in v] + elif isinstance(v, dict): + v = {key: str(val) for key, val in v.items()} else: v = str(v) + if isinstance(v, SpatialReferences): v = ' '.join(str(s) for s in v.references) or None + if isinstance(v, Reference): v = str(v) or None + out[k] = v return out @@ -533,7 +542,7 @@ def _process_value(value): class workflow(_Config): """Configure the particular execution graph of this workflow.""" - err_on_warn = None + err_on_warn = False """Cast AROMA warnings to errors.""" melodic_dim = None """Number of ICA components to be estimated by MELODIC @@ -729,8 +738,8 @@ def init_spaces(checkpoint=True): spaces.checkpoint() # Add the default standard space if not already present (required by several sub-workflows) - if 'MNI152NLin2009cAsym' not in spaces.get_spaces(nonstandard=False, dim=(3,)): - spaces.add(Reference('MNI152NLin2009cAsym', {})) + if 'MNI152NLin6Asym' not in spaces.get_spaces(nonstandard=False, dim=(3,)): + spaces.add(Reference('MNI152NLin6Asym', {})) # Ensure user-defined spatial references for outputs are correctly parsed. # Certain options require normalization to a space not explicitly defined by users. diff --git a/src/fmripost_aroma/data/boilerplate.bib b/src/fmripost_aroma/data/boilerplate.bib new file mode 100644 index 0000000..bbb15c5 --- /dev/null +++ b/src/fmripost_aroma/data/boilerplate.bib @@ -0,0 +1,388 @@ +@article{fmriprep1, + author = {Esteban, Oscar and Markiewicz, Christopher and Blair, Ross W and Moodie, Craig and Isik, Ayse Ilkay and Erramuzpe Aliaga, Asier and Kent, James and Goncalves, Mathias and DuPre, Elizabeth and Snyder, Madeleine and Oya, Hiroyuki and Ghosh, Satrajit and Wright, Jessey and Durnez, Joke and Poldrack, Russell and Gorgolewski, Krzysztof Jacek}, + title = {{fMRIPrep}: a robust preprocessing pipeline for functional {MRI}}, + volume = {16}, + pages = {111--116}, + year = {2019}, + doi = {10.1038/s41592-018-0235-4}, + journal = {Nature Methods} +} + +@article{fmriprep2, + author = {Esteban, Oscar and Blair, Ross and Markiewicz, Christopher J. and Berleant, Shoshana L. and Moodie, Craig and Ma, Feilong and Isik, Ayse Ilkay and Erramuzpe, Asier and Kent, James D. andGoncalves, Mathias and DuPre, Elizabeth and Sitek, Kevin R. and Gomez, Daniel E. P. and Lurie, Daniel J. and Ye, Zhifang and Poldrack, Russell A. and Gorgolewski, Krzysztof J.}, + title = {fMRIPrep }, + year = 2018, + doi = {10.5281/zenodo.852659}, + publisher = {Zenodo}, + journal = {Software} +} + +@article{templateflow, + author = {Ciric, R. and Thompson, William H. and Lorenz, R. and Goncalves, M. and MacNicol, E. and Markiewicz, C. J. and Halchenko, Y. O. and Ghosh, S. S. and Gorgolewski, K. J. and Poldrack, R. A. and Esteban, O.}, + title = {{TemplateFlow}: {FAIR}-sharing of multi-scale, multi-species brain models}, + volume = {19}, + pages = {1568--1571}, + year = {2022}, + doi = {10.1038/s41592-022-01681-2}, + journal = {Nature Methods} +} + +@article{nipype1, + author = {Gorgolewski, K. and Burns, C. D. and Madison, C. and Clark, D. and Halchenko, Y. O. and Waskom, M. L. and Ghosh, S.}, + doi = {10.3389/fninf.2011.00013}, + journal = {Frontiers in Neuroinformatics}, + pages = 13, + shorttitle = {Nipype}, + title = {Nipype: a flexible, lightweight and extensible neuroimaging data processing framework in Python}, + volume = 5, + year = 2011 +} + +@article{nipype2, + author = {Gorgolewski, Krzysztof J. and Esteban, Oscar and Markiewicz, Christopher J. and Ziegler, Erik and Ellis, David Gage and Notter, Michael Philipp and Jarecka, Dorota and Johnson, Hans and Burns, Christopher and Manhães-Savio, Alexandre and Hamalainen, Carlo and Yvernault, Benjamin and Salo, Taylor and Jordan, Kesshi and Goncalves, Mathias and Waskom, Michael and Clark, Daniel and Wong, Jason and Loney, Fred and Modat, Marc and Dewey, Blake E and Madison, Cindee and Visconti di Oleggio Castello, Matteo and Clark, Michael G. and Dayan, Michael and Clark, Dav and Keshavan, Anisha and Pinsard, Basile and Gramfort, Alexandre and Berleant, Shoshana and Nielson, Dylan M. and Bougacha, Salma and Varoquaux, Gael and Cipollini, Ben and Markello, Ross and Rokem, Ariel and Moloney, Brendan and Halchenko, Yaroslav O. and Wassermann , Demian and Hanke, Michael and Horea, Christian and Kaczmarzyk, Jakub and Gilles de Hollander and DuPre, Elizabeth and Gillman, Ashley and Mordom, David and Buchanan, Colin and Tungaraza, Rosalia and Pauli, Wolfgang M. and Iqbal, Shariq and Sikka, Sharad and Mancini, Matteo and Schwartz, Yannick and Malone, Ian B. and Dubois, Mathieu and Frohlich, Caroline and Welch, David and Forbes, Jessica and Kent, James and Watanabe, Aimi and Cumba, Chad and Huntenburg, Julia M. and Kastman, Erik and Nichols, B. Nolan and Eshaghi, Arman and Ginsburg, Daniel and Schaefer, Alexander and Acland, Benjamin and Giavasis, Steven and Kleesiek, Jens and Erickson, Drew and Küttner, René and Haselgrove, Christian and Correa, Carlos and Ghayoor, Ali and Liem, Franz and Millman, Jarrod and Haehn, Daniel and Lai, Jeff and Zhou, Dale and Blair, Ross and Glatard, Tristan and Renfro, Mandy and Liu, Siqi and Kahn, Ari E. and Pérez-García, Fernando and Triplett, William and Lampe, Leonie and Stadler, Jörg and Kong, Xiang-Zhen and Hallquist, Michael and Chetverikov, Andrey and Salvatore, John and Park, Anne and Poldrack, Russell and Craddock, R. Cameron and Inati, Souheil and Hinds, Oliver and Cooper, Gavin and Perkins, L. Nathan and Marina, Ana and Mattfeld, Aaron and Noel, Maxime and Lukas Snoek and Matsubara, K and Cheung, Brian and Rothmei, Simon and Urchs, Sebastian and Durnez, Joke and Mertz, Fred and Geisler, Daniel and Floren, Andrew and Gerhard, Stephan and Sharp, Paul and Molina-Romero, Miguel and Weinstein, Alejandro and Broderick, William and Saase, Victor and Andberg, Sami Kristian and Harms, Robbert and Schlamp, Kai and Arias, Jaime and Papadopoulos Orfanos, Dimitri and Tarbert, Claire and Tambini, Arielle and De La Vega, Alejandro and Nickson, Thomas and Brett, Matthew and Falkiewicz, Marcel and Podranski, Kornelius and Linkersdörfer, Janosch and Flandin, Guillaume and Ort, Eduard and Shachnev, Dmitry and McNamee, Daniel and Davison, Andrew and Varada, Jan and Schwabacher, Isaac and Pellman, John and Perez-Guevara, Martin and Khanuja, Ranjeet and Pannetier, Nicolas and McDermottroe, Conor and Ghosh, Satrajit}, + title = {Nipype}, + year = 2018, + doi = {10.5281/zenodo.596855}, + publisher = {Zenodo}, + journal = {Software} +} + +@article{n4, + author = {Tustison, N. J. and Avants, B. B. and Cook, P. A. and Zheng, Y. and Egan, A. and Yushkevich, P. A. and Gee, J. C.}, + doi = {10.1109/TMI.2010.2046908}, + issn = {0278-0062}, + journal = {IEEE Transactions on Medical Imaging}, + number = 6, + pages = {1310-1320}, + shorttitle = {N4ITK}, + title = {N4ITK: Improved N3 Bias Correction}, + volume = 29, + year = 2010 +} + +@article{fs_reconall, + author = {Dale, Anders M. and Fischl, Bruce and Sereno, Martin I.}, + doi = {10.1006/nimg.1998.0395}, + issn = {1053-8119}, + journal = {NeuroImage}, + number = 2, + pages = {179-194}, + shorttitle = {Cortical Surface-Based Analysis}, + title = {Cortical Surface-Based Analysis: I. Segmentation and Surface Reconstruction}, + url = {https://www.sciencedirect.com/science/article/pii/S1053811998903950}, + volume = 9, + year = 1999 +} + + + +@article{mindboggle, + author = {Klein, Arno and Ghosh, Satrajit S. and Bao, Forrest S. and Giard, Joachim and Häme, Yrjö and Stavsky, Eliezer and Lee, Noah and Rossa, Brian and Reuter, Martin and Neto, Elias Chaibub and Keshavan, Anisha}, + doi = {10.1371/journal.pcbi.1005350}, + issn = {1553-7358}, + journal = {PLOS Computational Biology}, + number = 2, + pages = {e1005350}, + title = {Mindboggling morphometry of human brains}, + url = {https://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1005350}, + volume = 13, + year = 2017 +} + +@article{mni152lin, + title = {A {Probabilistic} {Atlas} of the {Human} {Brain}: {Theory} and {Rationale} for {Its} {Development}: {The} {International} {Consortium} for {Brain} {Mapping} ({ICBM})}, + author = {Mazziotta, John C. and Toga, Arthur W. and Evans, Alan and Fox, Peter and Lancaster, Jack}, + volume = {2}, + issn = {1053-8119}, + shorttitle = {A {Probabilistic} {Atlas} of the {Human} {Brain}}, + doi = {10.1006/nimg.1995.1012}, + number = {2, Part A}, + journal = {NeuroImage}, + year = {1995}, + pages = {89--101} +} + +@article{mni152nlin2009casym, + title = {Unbiased nonlinear average age-appropriate brain templates from birth to adulthood}, + author = {Fonov, VS and Evans, AC and McKinstry, RC and Almli, CR and Collins, DL}, + doi = {10.1016/S1053-8119(09)70884-5}, + journal = {NeuroImage}, + pages = {S102}, + volume = {47, Supplement 1}, + year = 2009 +} + +@article{mni152nlin6asym, + author = {Evans, AC and Janke, AL and Collins, DL and Baillet, S}, + title = {Brain templates and atlases}, + doi = {10.1016/j.neuroimage.2012.01.024}, + journal = {NeuroImage}, + volume = {62}, + number = {2}, + pages = {911--922}, + year = 2012 +} + +@article{ants, + author = {Avants, B.B. and Epstein, C.L. and Grossman, M. and Gee, J.C.}, + doi = {10.1016/j.media.2007.06.004}, + issn = {1361-8415}, + journal = {Medical Image Analysis}, + number = 1, + pages = {26-41}, + shorttitle = {Symmetric diffeomorphic image registration with cross-correlation}, + title = {Symmetric diffeomorphic image registration with cross-correlation: Evaluating automated labeling of elderly and neurodegenerative brain}, + url = {https://www.sciencedirect.com/science/article/pii/S1361841507000606}, + volume = 12, + year = 2008 +} + +@article{fsl_fast, + author = {Zhang, Y. and Brady, M. and Smith, S.}, + doi = {10.1109/42.906424}, + issn = {0278-0062}, + journal = {IEEE Transactions on Medical Imaging}, + number = 1, + pages = {45-57}, + title = {Segmentation of brain {MR} images through a hidden Markov random field model and the expectation-maximization algorithm}, + volume = 20, + year = 2001 +} + + +@article{fieldmapless1, + author = {Wang, Sijia and Peterson, Daniel J. and Gatenby, J. C. and Li, Wenbin and Grabowski, Thomas J. and Madhyastha, Tara M.}, + doi = {10.3389/fninf.2017.00017}, + issn = {1662-5196}, + journal = {Frontiers in Neuroinformatics}, + language = {English}, + title = {Evaluation of Field Map and Nonlinear Registration Methods for Correction of Susceptibility Artifacts in Diffusion {MRI}}, + url = {https://www.frontiersin.org/articles/10.3389/fninf.2017.00017/full}, + volume = 11, + year = 2017 +} + +@phdthesis{fieldmapless2, + address = {Berlin}, + author = {Huntenburg, Julia M.}, + language = {eng}, + school = {Freie Universität}, + title = {Evaluating nonlinear coregistration of {BOLD} {EPI} and T1w images}, + type = {Master's Thesis}, + url = {https://hdl.handle.net/11858/00-001M-0000-002B-1CB5-A}, + year = 2014 +} + +@article{fieldmapless3, + author = {Treiber, Jeffrey Mark and White, Nathan S. and Steed, Tyler Christian and Bartsch, Hauke and Holland, Dominic and Farid, Nikdokht and McDonald, Carrie R. and Carter, Bob S. and Dale, Anders Martin and Chen, Clark C.}, + doi = {10.1371/journal.pone.0152472}, + issn = {1932-6203}, + journal = {PLOS ONE}, + number = 3, + pages = {e0152472}, + title = {Characterization and Correction of Geometric Distortions in 814 Diffusion Weighted Images}, + url = {https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0152472}, + volume = 11, + year = 2016 +} + +@article{flirt, + title = {A global optimisation method for robust affine registration of brain images}, + volume = {5}, + issn = {1361-8415}, + url = {https://www.sciencedirect.com/science/article/pii/S1361841501000366}, + doi = {10.1016/S1361-8415(01)00036-6}, + number = {2}, + urldate = {2018-07-27}, + journal = {Medical Image Analysis}, + author = {Jenkinson, Mark and Smith, Stephen}, + year = {2001}, + keywords = {Affine transformation, flirt, fsl, Global optimisation, Multi-resolution search, Multimodal registration, Robustness}, + pages = {143--156} +} + +@article{mcflirt, + author = {Jenkinson, Mark and Bannister, Peter and Brady, Michael and Smith, Stephen}, + doi = {10.1006/nimg.2002.1132}, + issn = {1053-8119}, + journal = {NeuroImage}, + number = 2, + pages = {825-841}, + title = {Improved Optimization for the Robust and Accurate Linear Registration and Motion Correction of Brain Images}, + url = {https://www.sciencedirect.com/science/article/pii/S1053811902911328}, + volume = 17, + year = 2002 +} + +@article{bbr, + author = {Greve, Douglas N and Fischl, Bruce}, + doi = {10.1016/j.neuroimage.2009.06.060}, + issn = {1095-9572}, + journal = {NeuroImage}, + number = 1, + pages = {63-72}, + title = {Accurate and robust brain image alignment using boundary-based registration}, + volume = 48, + year = 2009 +} + +@article{power_fd_dvars, + author = {Power, Jonathan D. and Mitra, Anish and Laumann, Timothy O. and Snyder, Abraham Z. and Schlaggar, Bradley L. and Petersen, Steven E.}, + doi = {10.1016/j.neuroimage.2013.08.048}, + issn = {1053-8119}, + journal = {NeuroImage}, + number = {Supplement C}, + pages = {320-341}, + title = {Methods to detect, characterize, and remove motion artifact in resting state fMRI}, + url = {https://www.sciencedirect.com/science/article/pii/S1053811913009117}, + volume = 84, + year = 2014 +} + +@article{confounds_satterthwaite_2013, + author = {Satterthwaite, Theodore D. and Elliott, Mark A. and Gerraty, Raphael T. and Ruparel, Kosha and Loughead, James and Calkins, Monica E. and Eickhoff, Simon B. and Hakonarson, Hakon and Gur, Ruben C. and Gur, Raquel E. and Wolf, Daniel H.}, + doi = {10.1016/j.neuroimage.2012.08.052}, + issn = {10538119}, + journal = {NeuroImage}, + number = 1, + pages = {240--256}, + title = {{An improved framework for confound regression and filtering for control of motion artifact in the preprocessing of resting-state functional connectivity data}}, + url = {https://www.sciencedirect.com/science/article/pii/S1053811912008609}, + volume = 64, + year = 2013 +} + + +@article{nilearn, + author = {Abraham, Alexandre and Pedregosa, Fabian and Eickenberg, Michael and Gervais, Philippe and Mueller, Andreas and Kossaifi, Jean and Gramfort, Alexandre and Thirion, Bertrand and Varoquaux, Gael}, + doi = {10.3389/fninf.2014.00014}, + issn = {1662-5196}, + journal = {Frontiers in Neuroinformatics}, + language = {English}, + title = {Machine learning for neuroimaging with scikit-learn}, + url = {https://www.frontiersin.org/articles/10.3389/fninf.2014.00014/full}, + volume = 8, + year = 2014 +} + +@article{lanczos, + author = {Lanczos, C.}, + doi = {10.1137/0701007}, + issn = {0887-459X}, + journal = {Journal of the Society for Industrial and Applied Mathematics Series B Numerical Analysis}, + number = 1, + pages = {76-85}, + title = {Evaluation of Noisy Data}, + url = {https://epubs.siam.org/doi/10.1137/0701007}, + volume = 1, + year = 1964 +} + +@article{compcor, + author = {Behzadi, Yashar and Restom, Khaled and Liau, Joy and Liu, Thomas T.}, + doi = {10.1016/j.neuroimage.2007.04.042}, + issn = {1053-8119}, + journal = {NeuroImage}, + number = 1, + pages = {90-101}, + title = {A component based noise correction method ({CompCor}) for {BOLD} and perfusion based fMRI}, + url = {https://www.sciencedirect.com/science/article/pii/S1053811907003837}, + volume = 37, + year = 2007 +} + +@article{hcppipelines, + author = {Glasser, Matthew F. and Sotiropoulos, Stamatios N. and Wilson, J. Anthony and Coalson, Timothy S. and Fischl, Bruce and Andersson, Jesper L. and Xu, Junqian and Jbabdi, Saad and Webster, Matthew and Polimeni, Jonathan R. and Van Essen, David C. and Jenkinson, Mark}, + doi = {10.1016/j.neuroimage.2013.04.127}, + issn = {1053-8119}, + journal = {NeuroImage}, + pages = {105-124}, + series = {Mapping the Connectome}, + title = {The minimal preprocessing pipelines for the Human Connectome Project}, + url = {https://www.sciencedirect.com/science/article/pii/S1053811913005053}, + volume = 80, + year = 2013 +} + +@article{fs_template, + author = {Reuter, Martin and Rosas, Herminia Diana and Fischl, Bruce}, + doi = {10.1016/j.neuroimage.2010.07.020}, + journal = {NeuroImage}, + number = 4, + pages = {1181-1196}, + title = {Highly accurate inverse consistent registration: A robust approach}, + volume = 53, + year = 2010 +} + +@article{afni, + author = {Cox, Robert W. and Hyde, James S.}, + doi = {10.1002/(SICI)1099-1492(199706/08)10:4/5<171::AID-NBM453>3.0.CO;2-L}, + journal = {NMR in Biomedicine}, + number = {4-5}, + pages = {171-178}, + title = {Software tools for analysis and visualization of fMRI data}, + volume = 10, + year = 1997 +} + +@article{posse_t2s, + author = {Posse, Stefan and Wiese, Stefan and Gembris, Daniel and Mathiak, Klaus and Kessler, Christoph and Grosse-Ruyken, Maria-Liisa and Elghahwagi, Barbara and Richards, Todd and Dager, Stephen R. and Kiselev, Valerij G.}, + doi = {10.1002/(SICI)1522-2594(199907)42:1<87::AID-MRM13>3.0.CO;2-O}, + journal = {Magnetic Resonance in Medicine}, + number = 1, + pages = {87-97}, + title = {Enhancement of {BOLD}-contrast sensitivity by single-shot multi-echo functional {MR} imaging}, + volume = 42, + year = 1999 +} + +@article{topup, + author = {Jesper L.R. Andersson and Stefan Skare and John Ashburner}, + title = {How to correct susceptibility distortions in spin-echo echo-planar images: application to diffusion tensor imaging}, + journal = {NeuroImage}, + volume = 20, + number = 2, + pages = {870-888}, + year = 2003, + issn = {1053-8119}, + doi = {10.1016/S1053-8119(03)00336-7}, + url = {https://www.sciencedirect.com/science/article/pii/S1053811903003367} +} + +@article{patriat_improved_2017, + title = {An improved model of motion-related signal changes in {fMRI}}, + volume = {144, Part A}, + issn = {1053-8119}, + url = {https://www.sciencedirect.com/science/article/pii/S1053811916304360}, + doi = {10.1016/j.neuroimage.2016.08.051}, + abstract = {Head motion is a significant source of noise in the estimation of functional connectivity from resting-state functional MRI (rs-fMRI). Current strategies to reduce this noise include image realignment, censoring time points corrupted by motion, and including motion realignment parameters and their derivatives as additional nuisance regressors in the general linear model. However, this nuisance regression approach assumes that the motion-induced signal changes are linearly related to the estimated realignment parameters, which is not always the case. In this study we develop an improved model of motion-related signal changes, where nuisance regressors are formed by first rotating and translating a single brain volume according to the estimated motion, re-registering the data, and then performing a principal components analysis (PCA) on the resultant time series of both moved and re-registered data. We show that these “Motion Simulated (MotSim)” regressors account for significantly greater fraction of variance, result in higher temporal signal-to-noise, and lead to functional connectivity estimates that are less affected by motion compared to the most common current approach of using the realignment parameters and their derivatives as nuisance regressors. This improvement should lead to more accurate estimates of functional connectivity, particularly in populations where motion is prevalent, such as patients and young children.}, + urldate = {2017-01-18}, + journal = {NeuroImage}, + author = {Patriat, Rémi and Reynolds, Richard C. and Birn, Rasmus M.}, + month = jan, + year = {2017}, + keywords = {Motion, Correction, Methods, Rs-fMRI}, + pages = {74--82}, +} + +@article{ica_aroma, + title={ICA-AROMA: A robust ICA-based strategy for removing motion artifacts from fMRI data}, + author={Pruim, Raimon HR and Mennes, Maarten and van Rooij, Daan and Llera, Alberto and Buitelaar, Jan K and Beckmann, Christian F}, + journal={Neuroimage}, + volume={112}, + pages={267--277}, + year={2015}, + publisher={Elsevier} +} + +@article{susan, + title={SUSAN—a new approach to low level image processing}, + author={Smith, Stephen M and Brady, J Michael}, + journal={International journal of computer vision}, + volume={23}, + number={1}, + pages={45--78}, + year={1997}, + publisher={Springer} +} diff --git a/src/fmripost_aroma/data/io_spec.json b/src/fmripost_aroma/data/io_spec.json index 400ea31..cb1b2e4 100644 --- a/src/fmripost_aroma/data/io_spec.json +++ b/src/fmripost_aroma/data/io_spec.json @@ -3,8 +3,11 @@ "raw": { "bold_raw": { "datatype": "func", - "space": null, - "desc": null, + "echo": null, + "part": [ + "mag", + null + ], "suffix": "bold", "extension": [ ".nii.gz", @@ -13,9 +16,14 @@ } }, "derivatives": { - "bold_boldref": { + "bold_std": { "datatype": "func", - "space": null, + "echo": null, + "part": [ + "mag", + null + ], + "space": "MNI152NLin6Asym", "desc": "preproc", "suffix": "bold", "extension": [ @@ -23,15 +31,51 @@ ".nii" ] }, - "bold_MNI152NLin6Asym": { + "bold_mask_std": { "datatype": "func", + "echo": null, + "part": [ + "mag", + null + ], "space": "MNI152NLin6Asym", - "desc": "preproc", - "suffix": "bold", + "desc": "brain", + "suffix": "mask", + "extension": [ + ".nii.gz", + ".nii" + ] + }, + "bold_mask": { + "datatype": "func", + "echo": null, + "part": [ + "mag", + null + ], + "space": null, + "desc": "brain", + "suffix": "mask", "extension": [ ".nii.gz", ".nii" ] + }, + "confounds": { + "datatype": "func", + "echo": null, + "part": [ + "mag", + null + ], + "space": null, + "res": null, + "den": null, + "desc": "confounds", + "suffix": "timeseries", + "extension": [ + ".tsv" + ] } }, "transforms": { @@ -46,7 +90,7 @@ "boldref2anat": { "datatype": "func", "from": "boldref", - "to": "anat", + "to": ["anat", "T1w", "T2w"], "mode": "image", "suffix": "xfm", "extension": ".txt" @@ -60,11 +104,13 @@ }, "anat2mni152nlin6asym": { "datatype": "anat", - "from": "anat", + "task": null, + "run": null, + "from": ["anat", "T1w", "T2w"], "to": "MNI152NLin6Asym", "mode": "image", "suffix": "xfm", - "extension": ".txt" + "extension": ".h5" } } }, @@ -92,13 +138,13 @@ } ], "patterns": [ - "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_res-{res}][_echo-{echo}][_part-{part}][_space-{space}][_desc-{desc}]_{suffix}.{extension|nii.gz}", - "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_res-{res}][_echo-{echo}][_part-{part}][_space-{space}][_stat-{statistic}][_desc-{desc}]_{suffix}.{extension|nii.gz}", - "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_res-{res}][_echo-{echo}][_part-{part}][_space-{space}][_stat-{statistic}][_desc-{desc}]_{suffix}.{extension|tsv}", + "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_res-{resolution}][_echo-{echo}][_part-{part}][_space-{space}][_desc-{desc}]_{suffix}.{extension|nii.gz}", + "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_res-{resolution}][_echo-{echo}][_part-{part}][_space-{space}][_stat-{statistic}][_desc-{desc}]_{suffix}.{extension|nii.gz}", + "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_res-{resolution}][_echo-{echo}][_part-{part}][_space-{space}][_stat-{statistic}][_desc-{desc}]_{suffix}.{extension|tsv}", "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_part-{part}][_desc-{desc}]_{suffix}.{extension|tsv}", "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_run-{run}][_part-{part}][_desc-{desc}]_{suffix}.{extension}", - "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_space-{space}][_res-{res}][_den-{den}][_hemi-{hemi}][_label-{label}][_desc-{desc}]_{suffix<|boldref|dseg|mask>}.{extension}", + "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_space-{space}][_res-{resolution}][_den-{density}][_hemi-{hemi}][_label-{label}][_desc-{desc}]_{suffix<|boldref|dseg|mask>}.{extension}", "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_run-{run}]_from-{from}_to-{to}_mode-{mode|image}_{suffix|xfm}.{extension}", - "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_res-{res}][_echo-{echo}][_part-{part}][_space-{space}][_stat-{statistic}][_desc-{desc}]_{suffix}.{extension|svg}" + "sub-{subject}[/ses-{session}]/{datatype|func}/sub-{subject}[_ses-{session}]_task-{task}[_acq-{acquisition}][_ce-{ceagent}][_rec-{reconstruction}][_res-{resolution}][_echo-{echo}][_part-{part}][_space-{space}][_stat-{statistic}][_desc-{desc}]_{suffix}.{extension|svg}" ] } diff --git a/src/fmripost_aroma/data/reports-spec.yml b/src/fmripost_aroma/data/reports-spec.yml new file mode 100644 index 0000000..a06f803 --- /dev/null +++ b/src/fmripost_aroma/data/reports-spec.yml @@ -0,0 +1,47 @@ +package: fmripost-aroma +sections: +- name: Summary + reportlets: + - bids: {datatype: figures, desc: summary, suffix: bold} + +- name: Functional + ordering: session,task,acquisition,ceagent,reconstruction,direction,run,echo + reportlets: + - bids: {datatype: figures, desc: summary, suffix: bold} + - bids: {datatype: figures, desc: validation, suffix: bold} + - bids: {datatype: figures, desc: aroma, suffix: bold} + - bids: {datatype: figures, desc: coreg, suffix: bold} + caption: This panel shows the alignment of the reference EPI (BOLD) image to the + anatomical (T1-weighted) image. + The reference EPI has been contrast enhanced and susceptibility-distortion + corrected (if applicable) for improved anatomical fidelity. + The anatomical image has been resampled into EPI space, as well as the + anatomical white matter mask, which appears as a red contour. + static: false + subtitle: Alignment of functional and anatomical MRI data (coregistration) + - bids: + datatype: figures + desc: '[aggr,nonaggr,orthaggr]carpetplot' + extension: [.html] + suffix: bold +- name: About + nested: true + reportlets: + - bids: {datatype: figures, desc: about, suffix: T1w} + - custom: boilerplate + path: '{out_dir}/logs' + bibfile: ['fmripost_aroma', 'data/boilerplate.bib'] + caption: | +

We kindly ask to report results preprocessed with this tool using the following boilerplate.

+ + title: Methods + - custom: errors + path: '{out_dir}/sub-{subject}/log/{run_uuid}' + captions: NiReports may have recorded failure conditions. + title: Errors diff --git a/src/fmripost_aroma/data/tests/config.toml b/src/fmripost_aroma/data/tests/config.toml index c08ed92..6f92048 100644 --- a/src/fmripost_aroma/data/tests/config.toml +++ b/src/fmripost_aroma/data/tests/config.toml @@ -13,8 +13,6 @@ aggr_ses_reports = false bids_dir = "/home/runner/work/fmripost-aroma/fmripost-aroma/src/fmripost_aroma/tests/data" bids_description_hash = "5d42e27751bbc884eca87cb4e62b9a0cca0cd86f8e578747fe89b77e6c5b21e5" boilerplate_only = false -fs_license_file = "/opt/freesurfer/license.txt" -fs_subjects_dir = "/opt/freesurfer/subjects" log_dir = "/home/runner/work/fmripost-aroma/fmripost-aroma/src/fmripost_aroma/tests/data" log_level = 40 low_mem = false @@ -22,7 +20,7 @@ md_only_boilerplate = false notrack = true fmripost_aroma_dir = "/home/runner/work/fmripost-aroma/fmripost-aroma/src/fmripost_aroma/tests/out" output_dir = "/home/runner/work/fmripost-aroma/fmripost-aroma/src/fmripost_aroma/tests/out" -output_spaces = "MNI152NLin2009cAsym:res-2 MNI152NLin2009cAsym:res-native fsaverage:den-10k fsaverage:den-30k" +output_spaces = "MNI152NLin6Asym:res-2" reports_only = false run_uuid = "20200306-105302_d365772b-fd60-4741-a722-372c2f558b50" participant_label = [ "01",] diff --git a/src/fmripost_aroma/interfaces/confounds.py b/src/fmripost_aroma/interfaces/confounds.py index 63fc168..580aab7 100644 --- a/src/fmripost_aroma/interfaces/confounds.py +++ b/src/fmripost_aroma/interfaces/confounds.py @@ -109,7 +109,14 @@ class _ICADenoiseInputSpec(BaseInterfaceInputSpec): method = traits.Enum('aggr', 'nonaggr', 'orthaggr', mandatory=True, desc='denoising method') bold_file = File(exists=True, mandatory=True, desc='input file to denoise') mask_file = File(exists=True, mandatory=True, desc='mask file') - confounds = File(exists=True, mandatory=True, desc='confounds file') + mixing = File(exists=True, mandatory=True, desc='mixing matrix file') + metrics = File(exists=True, mandatory=True, desc='metrics file') + skip_vols = traits.Int( + desc=( + 'Number of non steady state volumes identified. ' + 'Will be removed from mixing, but not bold' + ), + ) class _ICADenoiseOutputSpec(TraitedSpec): @@ -130,17 +137,19 @@ def _run_interface(self, runtime): method = self.inputs.method bold_file = self.inputs.bold_file - confounds_file = self.inputs.confounds - metrics_file = self.inputs.metrics_file + mixing_file = self.inputs.mixing + metrics_file = self.inputs.metrics - confounds_df = pd.read_table(confounds_file) + # Load the mixing matrix. It doesn't have a header row + mixing = np.loadtxt(mixing_file) + mixing = mixing[self.inputs.skip_vols :, :] # Split up component time series into accepted and rejected components metrics_df = pd.read_table(metrics_file) - rejected_columns = metrics_df.loc[metrics_df['classification'] == 'rejected', 'Component'] - accepted_columns = metrics_df.loc[metrics_df['classification'] == 'accepted', 'Component'] - rejected_components = confounds_df[rejected_columns].to_numpy() - accepted_components = confounds_df[accepted_columns].to_numpy() + rejected_idx = metrics_df.loc[metrics_df['classification'] == 'rejected'].index.values + accepted_idx = metrics_df.loc[metrics_df['classification'] == 'accepted'].index.values + rejected_components = mixing[:, rejected_idx] + accepted_components = mixing[:, accepted_idx] # Z-score all of the components rejected_components = stats.zscore(rejected_components, axis=0) accepted_components = stats.zscore(accepted_components, axis=0) @@ -198,7 +207,7 @@ def _run_interface(self, runtime): ( rejected_components, accepted_components, - np.ones((confounds_df.shape[0], 1)), + np.ones((mixing.shape[0], 1)), ), ) betas = np.linalg.lstsq(regressors, data, rcond=None)[0][:-1] diff --git a/src/fmripost_aroma/reports/core.py b/src/fmripost_aroma/reports/core.py new file mode 100644 index 0000000..4bf2ab2 --- /dev/null +++ b/src/fmripost_aroma/reports/core.py @@ -0,0 +1,148 @@ +# emacs: -*- mode: python; py-indent-offset: 4; indent-tabs-mode: nil -*- +# vi: set ft=python sts=4 ts=4 sw=4 et: +# +# Copyright The NiPreps Developers +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. +# +# We support and encourage derived works from this project, please read +# about our expectations at +# +# https://www.nipreps.org/community/licensing/ +# +from pathlib import Path + +from fmripost_aroma import config, data +from nireports.assembler.report import Report + + +def run_reports( + output_dir, + subject_label, + run_uuid, + bootstrap_file=None, + out_filename='report.html', + reportlets_dir=None, + errorname='report.err', + **entities, +): + """Run the reports.""" + robj = Report( + output_dir, + run_uuid, + bootstrap_file=bootstrap_file, + out_filename=out_filename, + reportlets_dir=reportlets_dir, + plugins=None, + plugin_meta=None, + metadata=None, + **entities, + ) + + # Count nbr of subject for which report generation failed + try: + robj.generate_report() + except: # noqa: E722 + import sys + import traceback + + # Store the list of subjects for which report generation failed + traceback.print_exception(*sys.exc_info(), file=str(Path(output_dir) / 'logs' / errorname)) + return subject_label + + return None + + +def generate_reports( + subject_list, + output_dir, + run_uuid, + session_list=None, + bootstrap_file=None, + work_dir=None, +): + """Generate reports for a list of subjects.""" + reportlets_dir = None + if work_dir is not None: + reportlets_dir = Path(work_dir) / 'reportlets' + + if isinstance(subject_list, str): + subject_list = [subject_list] + + errors = [] + for subject_label in subject_list: + # The number of sessions is intentionally not based on session_list but + # on the total number of sessions, because I want the final derivatives + # folder to be the same whether sessions were run one at a time or all-together. + n_ses = len(config.execution.layout.get_sessions(subject=subject_label)) + + if bootstrap_file is not None: + # If a config file is precised, we do not override it + html_report = 'report.html' + elif n_ses <= config.execution.aggr_ses_reports: + # If there are only a few session for this subject, + # we aggregate them in a single visual report. + bootstrap_file = data.load('reports-spec.yml') + html_report = 'report.html' + else: + # Beyond a threshold, we separate the anatomical report from the functional. + bootstrap_file = data.load('reports-spec-anat.yml') + html_report = f'sub-{subject_label.lstrip("sub-")}_anat.html' + + report_error = run_reports( + output_dir, + subject_label, + run_uuid, + bootstrap_file=bootstrap_file, + out_filename=html_report, + reportlets_dir=reportlets_dir, + errorname=f'report-{run_uuid}-{subject_label}.err', + subject=subject_label, + ) + # If the report generation failed, append the subject label for which it failed + if report_error is not None: + errors.append(report_error) + + if n_ses > config.execution.aggr_ses_reports: + # Beyond a certain number of sessions per subject, + # we separate the functional reports per session + if session_list is None: + all_filters = config.execution.bids_filters or {} + filters = all_filters.get('bold', {}) + session_list = config.execution.layout.get_sessions( + subject=subject_label, **filters + ) + + # Drop ses- prefixes + session_list = [ses[4:] if ses.startswith('ses-') else ses for ses in session_list] + + for session_label in session_list: + bootstrap_file = data.load('reports-spec-func.yml') + html_report = f'sub-{subject_label.lstrip("sub-")}_ses-{session_label}_func.html' + + report_error = run_reports( + output_dir, + subject_label, + run_uuid, + bootstrap_file=bootstrap_file, + out_filename=html_report, + reportlets_dir=reportlets_dir, + errorname=f'report-{run_uuid}-{subject_label}-func.err', + subject=subject_label, + session=session_label, + ) + # If the report generation failed, append the subject label for which it failed + if report_error is not None: + errors.append(report_error) + + return errors diff --git a/src/fmripost_aroma/tests/conftest.py b/src/fmripost_aroma/tests/conftest.py index 98738fe..f13243c 100644 --- a/src/fmripost_aroma/tests/conftest.py +++ b/src/fmripost_aroma/tests/conftest.py @@ -3,6 +3,7 @@ import base64 import importlib.resources import os +import re import pytest @@ -44,3 +45,209 @@ def base_parser(): parser = ArgumentParser(description='Test parser') return parser + + +@pytest.fixture(scope='session') +def base_ignore_list(): + """Create the standard ignore list used by fMRIPost-AROMA.""" + return [ + 'code', + 'stimuli', + 'sourcedata', + 'models', + re.compile(r'^\.'), + re.compile(r'sub-[a-zA-Z0-9]+(/ses-[a-zA-Z0-9]+)?/(beh|dwi|eeg|ieeg|meg|perf)'), + ] + + +@pytest.fixture(scope='session') +def minimal_ignore_list(base_ignore_list): + """Create an ignore list that ignores full derivative files from ds000005.""" + base_ignore_list = base_ignore_list[:] + files_to_ignore = [ + 'sub-01/anat/sub-01_desc-brain_mask.json', + 'sub-01/anat/sub-01_desc-brain_mask.nii.gz', + 'sub-01/anat/sub-01_desc-preproc_T1w.json', + 'sub-01/anat/sub-01_desc-preproc_T1w.nii.gz', + 'sub-01/anat/sub-01_desc-ribbon_mask.json', + 'sub-01/anat/sub-01_desc-ribbon_mask.nii.gz', + 'sub-01/anat/sub-01_dseg.nii.gz', + # 'sub-01/anat/sub-01_from-MNI152NLin2009cAsym_to-T1w_mode-image_xfm.h5', + # 'sub-01/anat/sub-01_from-MNI152NLin6Asym_to-T1w_mode-image_xfm.h5', + # 'sub-01/anat/sub-01_from-T1w_to-MNI152NLin2009cAsym_mode-image_xfm.h5', + # 'sub-01/anat/sub-01_from-T1w_to-MNI152NLin6Asym_mode-image_xfm.h5', + # 'sub-01/anat/sub-01_from-T1w_to-fsnative_mode-image_xfm.txt', + # 'sub-01/anat/sub-01_from-fsnative_to-T1w_mode-image_xfm.txt', + 'sub-01/anat/sub-01_hemi-L_desc-reg_sphere.surf.gii', + 'sub-01/anat/sub-01_hemi-L_midthickness.surf.gii', + 'sub-01/anat/sub-01_hemi-L_pial.surf.gii', + 'sub-01/anat/sub-01_hemi-L_space-fsLR_desc-msmsulc_sphere.surf.gii', + 'sub-01/anat/sub-01_hemi-L_space-fsLR_desc-reg_sphere.surf.gii', + 'sub-01/anat/sub-01_hemi-L_sphere.surf.gii', + 'sub-01/anat/sub-01_hemi-L_sulc.shape.gii', + 'sub-01/anat/sub-01_hemi-L_thickness.shape.gii', + 'sub-01/anat/sub-01_hemi-L_white.surf.gii', + 'sub-01/anat/sub-01_hemi-R_desc-reg_sphere.surf.gii', + 'sub-01/anat/sub-01_hemi-R_midthickness.surf.gii', + 'sub-01/anat/sub-01_hemi-R_pial.surf.gii', + 'sub-01/anat/sub-01_hemi-R_space-fsLR_desc-msmsulc_sphere.surf.gii', + 'sub-01/anat/sub-01_hemi-R_space-fsLR_desc-reg_sphere.surf.gii', + 'sub-01/anat/sub-01_hemi-R_sphere.surf.gii', + 'sub-01/anat/sub-01_hemi-R_sulc.shape.gii', + 'sub-01/anat/sub-01_hemi-R_thickness.shape.gii', + 'sub-01/anat/sub-01_hemi-R_white.surf.gii', + 'sub-01/anat/sub-01_label-CSF_probseg.nii.gz', + 'sub-01/anat/sub-01_label-GM_probseg.nii.gz', + 'sub-01/anat/sub-01_label-WM_probseg.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin2009cAsym_desc-brain_mask.json', + 'sub-01/anat/sub-01_space-MNI152NLin2009cAsym_desc-brain_mask.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin2009cAsym_desc-preproc_T1w.json', + 'sub-01/anat/sub-01_space-MNI152NLin2009cAsym_desc-preproc_T1w.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin2009cAsym_dseg.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin2009cAsym_label-CSF_probseg.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin2009cAsym_label-GM_probseg.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin2009cAsym_label-WM_probseg.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin6Asym_res-2_desc-brain_mask.json', + 'sub-01/anat/sub-01_space-MNI152NLin6Asym_res-2_desc-brain_mask.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin6Asym_res-2_desc-preproc_T1w.json', + 'sub-01/anat/sub-01_space-MNI152NLin6Asym_res-2_desc-preproc_T1w.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin6Asym_res-2_dseg.json', + 'sub-01/anat/sub-01_space-MNI152NLin6Asym_res-2_dseg.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin6Asym_res-2_label-CSF_probseg.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin6Asym_res-2_label-GM_probseg.nii.gz', + 'sub-01/anat/sub-01_space-MNI152NLin6Asym_res-2_label-WM_probseg.nii.gz', + 'sub-01/figures/', + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_desc-confounds_timeseries.json', + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_desc-confounds_timeseries.tsv', + # 'sub-01/func/sub-01_task-mixedgamblestask_run-01_desc-coreg_boldref.json', + # 'sub-01/func/sub-01_task-mixedgamblestask_run-01_desc-coreg_boldref.nii.gz', + # 'sub-01/func/sub-01_task-mixedgamblestask_run-01_desc-hmc_boldref.json', + # 'sub-01/func/sub-01_task-mixedgamblestask_run-01_desc-hmc_boldref.nii.gz', + # ( + # 'sub-01/func/sub-01_task-mixedgamblestask_run-01_from-boldref_to-T1w_mode-image_' + # 'desc-coreg_xfm.json' + # ), + # ( + # 'sub-01/func/sub-01_task-mixedgamblestask_run-01_from-boldref_to-T1w_mode-image_' + # 'desc-coreg_xfm.txt' + # ), + # ( + # 'sub-01/func/sub-01_task-mixedgamblestask_run-01_from-orig_to-boldref_mode-image_' + # 'desc-hmc_xfm.json' + # ), + # ( + # 'sub-01/func/sub-01_task-mixedgamblestask_run-01_from-orig_to-boldref_mode-image_' + # 'desc-hmc_xfm.txt' + # ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_space-MNI152NLin2009cAsym_' + 'boldref.nii.gz' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_space-MNI152NLin2009cAsym_' + 'desc-brain_mask.nii.gz' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_space-MNI152NLin2009cAsym_' + 'desc-preproc_bold.json' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_space-MNI152NLin2009cAsym_' + 'desc-preproc_bold.nii.gz' + ), + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_space-MNI152NLin6Asym_res-2_boldref.json', + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_space-MNI152NLin6Asym_' + 'res-2_boldref.nii.gz' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_space-MNI152NLin6Asym_' + 'res-2_desc-brain_mask.json' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_space-MNI152NLin6Asym_' + 'res-2_desc-brain_mask.nii.gz' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_space-MNI152NLin6Asym_' + 'res-2_desc-preproc_bold.json' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_space-MNI152NLin6Asym_' + 'res-2_desc-preproc_bold.nii.gz' + ), + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_desc-confounds_timeseries.json', + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_desc-confounds_timeseries.tsv', + # 'sub-01/func/sub-01_task-mixedgamblestask_run-02_desc-coreg_boldref.json', + # 'sub-01/func/sub-01_task-mixedgamblestask_run-02_desc-coreg_boldref.nii.gz', + # 'sub-01/func/sub-01_task-mixedgamblestask_run-02_desc-hmc_boldref.json', + # 'sub-01/func/sub-01_task-mixedgamblestask_run-02_desc-hmc_boldref.nii.gz', + # ( + # 'sub-01/func/sub-01_task-mixedgamblestask_run-02_from-boldref_to-T1w_mode-image_' + # 'desc-coreg_xfm.json' + # ), + # ( + # 'sub-01/func/sub-01_task-mixedgamblestask_run-02_from-boldref_to-T1w_mode-image_' + # 'desc-coreg_xfm.txt' + # ), + # ( + # 'sub-01/func/sub-01_task-mixedgamblestask_run-02_from-orig_to-boldref_mode-image_' + # 'desc-hmc_xfm.json' + # ), + # ( + # 'sub-01/func/sub-01_task-mixedgamblestask_run-02_from-orig_to-boldref_mode-image_' + # 'desc-hmc_xfm.txt' + # ), + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_space-MNI152NLin2009cAsym_boldref.nii.gz', + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_space-MNI152NLin2009cAsym_' + 'desc-brain_mask.nii.gz' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_space-MNI152NLin2009cAsym_' + 'desc-preproc_bold.json' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_space-MNI152NLin2009cAsym_' + 'desc-preproc_bold.nii.gz' + ), + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_space-MNI152NLin6Asym_res-2_boldref.json', + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_space-MNI152NLin6Asym_' + 'res-2_boldref.nii.gz' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_space-MNI152NLin6Asym_' + 'res-2_desc-brain_mask.json' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_space-MNI152NLin6Asym_' + 'res-2_desc-brain_mask.nii.gz' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_space-MNI152NLin6Asym_' + 'res-2_desc-preproc_bold.json' + ), + ( + 'sub-01/func/sub-01_task-mixedgamblestask_run-02_space-MNI152NLin6Asym_' + 'res-2_desc-preproc_bold.nii.gz' + ), + ] + regex = '|'.join(files_to_ignore) + base_ignore_list.append(re.compile(regex)) + return base_ignore_list + + +@pytest.fixture(scope='session') +def full_ignore_list(base_ignore_list): + """Create an ignore list that ignores minimal derivative files from ds000005. + + The 'full' run didn't include func-space outputs, so there's no func-space brain mask. + """ + base_ignore_list = base_ignore_list[:] + files_to_ignore = [ + 'sub-01/func/sub-01_task-mixedgamblestask_run-01_desc-brain_mask.nii.gz', + ] + regex = '|'.join(files_to_ignore) + base_ignore_list.append(re.compile(regex)) + return base_ignore_list diff --git a/src/fmripost_aroma/tests/test_base.py b/src/fmripost_aroma/tests/test_base.py index 9dc8e64..4f1b7bd 100644 --- a/src/fmripost_aroma/tests/test_base.py +++ b/src/fmripost_aroma/tests/test_base.py @@ -10,8 +10,8 @@ def test_init_ica_aroma_wf(tmp_path_factory): tempdir = tmp_path_factory.mktemp('test_init_ica_aroma_wf') with mock_config(): - config.workflow.fmripost_aroma_dir = tempdir / 'out' - config.workflow.work_dir = tempdir / 'work' + config.execution.fmripost_aroma_dir = tempdir / 'out' + config.execution.work_dir = tempdir / 'work' config.workflow.denoise_method = ['nonaggr', 'orthaggr'] config.workflow.melodic_dim = -200 config.workflow.err_on_warn = False @@ -29,8 +29,8 @@ def test_init_denoise_wf(tmp_path_factory): tempdir = tmp_path_factory.mktemp('test_init_denoise_wf') with mock_config(): - config.workflow.fmripost_aroma_dir = tempdir / 'out' - config.workflow.work_dir = tempdir / 'work' + config.execution.fmripost_aroma_dir = tempdir / 'out' + config.execution.work_dir = tempdir / 'work' wf = init_denoise_wf(bold_file='sub-01_task-rest_bold.nii.gz') assert wf.name == 'denoise_task_rest_wf' diff --git a/src/fmripost_aroma/tests/test_utils_bids.py b/src/fmripost_aroma/tests/test_utils_bids.py new file mode 100644 index 0000000..1742908 --- /dev/null +++ b/src/fmripost_aroma/tests/test_utils_bids.py @@ -0,0 +1,150 @@ +"""Lightweight tests for fmripost_aroma.utils.bids.""" + +import os + +import pytest +from bids.layout import BIDSLayout, BIDSLayoutIndexer +from fmripost_aroma.tests.utils import get_test_data_path +from fmripost_aroma.utils import bids as xbids + + +def test_collect_derivatives_raw(base_ignore_list): + """Test collect_derivatives with a raw dataset.""" + data_dir = get_test_data_path() + + raw_dataset = data_dir / 'ds000005-fmriprep' / 'sourcedata' / 'raw' + raw_layout = BIDSLayout( + raw_dataset, + config=['bids'], + validate=False, + indexer=BIDSLayoutIndexer(validate=False, index_metadata=False, ignore=base_ignore_list), + ) + + subject_data = xbids.collect_derivatives( + raw_dataset=raw_layout, + derivatives_dataset=None, + entities={'subject': '01', 'task': 'mixedgamblestask'}, + fieldmap_id=None, + spec=None, + patterns=None, + allow_multiple=True, + ) + expected = { + 'bold_raw': [ + 'sub-01_task-mixedgamblestask_run-01_bold.nii.gz', + 'sub-01_task-mixedgamblestask_run-02_bold.nii.gz', + ], + } + check_expected(subject_data, expected) + + with pytest.raises(ValueError, match='Multiple files found'): + xbids.collect_derivatives( + raw_dataset=raw_layout, + derivatives_dataset=None, + entities={'subject': '01', 'task': 'mixedgamblestask'}, + fieldmap_id=None, + spec=None, + patterns=None, + allow_multiple=False, + ) + + +def test_collect_derivatives_minimal(minimal_ignore_list): + """Test collect_derivatives with a minimal-mode dataset.""" + data_dir = get_test_data_path() + + derivatives_dataset = data_dir / 'ds000005-fmriprep' + derivatives_layout = BIDSLayout( + derivatives_dataset, + config=['bids', 'derivatives'], + validate=False, + indexer=BIDSLayoutIndexer( + validate=False, + index_metadata=False, + ignore=minimal_ignore_list, + ), + ) + + subject_data = xbids.collect_derivatives( + raw_dataset=None, + derivatives_dataset=derivatives_layout, + entities={'subject': '01', 'task': 'mixedgamblestask', 'run': 1}, + fieldmap_id=None, + spec=None, + patterns=None, + ) + expected = { + 'bold_std': None, + 'bold_mask_std': None, + # TODO: Add bold_mask to the dataset + # 'bold_mask': 'sub-01_task-mixedgamblestask_run-01_desc-brain_mask.nii.gz', + 'bold_mask': None, + 'confounds': None, + 'hmc': ( + 'sub-01_task-mixedgamblestask_run-01_from-orig_to-boldref_mode-image_desc-hmc_xfm.txt' + ), + 'boldref2anat': ( + 'sub-01_task-mixedgamblestask_run-01_from-boldref_to-T1w_mode-image_desc-coreg_xfm.txt' + ), + 'boldref2fmap': None, + 'anat2mni152nlin6asym': 'sub-01_from-T1w_to-MNI152NLin6Asym_mode-image_xfm.h5', + } + check_expected(subject_data, expected) + + +def test_collect_derivatives_full(full_ignore_list): + """Test collect_derivatives with a full-mode dataset.""" + data_dir = get_test_data_path() + + derivatives_dataset = data_dir / 'ds000005-fmriprep' + derivatives_layout = BIDSLayout( + derivatives_dataset, + config=['bids', 'derivatives'], + validate=False, + indexer=BIDSLayoutIndexer(validate=False, index_metadata=False, ignore=full_ignore_list), + ) + + subject_data = xbids.collect_derivatives( + raw_dataset=None, + derivatives_dataset=derivatives_layout, + entities={'subject': '01', 'task': 'mixedgamblestask', 'run': 1}, + fieldmap_id=None, + spec=None, + patterns=None, + ) + expected = { + 'bold_std': ( + 'sub-01_task-mixedgamblestask_run-01_space-MNI152NLin6Asym_res-2_desc-preproc_' + 'bold.nii.gz' + ), + 'bold_mask_std': ( + 'sub-01_task-mixedgamblestask_run-01_space-MNI152NLin6Asym_res-2_desc-brain_' + 'mask.nii.gz' + ), + 'bold_mask': None, + 'confounds': 'sub-01_task-mixedgamblestask_run-01_desc-confounds_timeseries.tsv', + 'hmc': ( + 'sub-01_task-mixedgamblestask_run-01_from-orig_to-boldref_mode-image_desc-hmc_xfm.txt' + ), + 'boldref2anat': ( + 'sub-01_task-mixedgamblestask_run-01_from-boldref_to-T1w_mode-image_desc-coreg_xfm.txt' + ), + 'boldref2fmap': None, + 'anat2mni152nlin6asym': 'sub-01_from-T1w_to-MNI152NLin6Asym_mode-image_xfm.h5', + } + check_expected(subject_data, expected) + + +def check_expected(subject_data, expected): + """Check expected values.""" + for key, value in expected.items(): + if isinstance(value, str): + assert subject_data[key] is not None, f'Key {key} is None.' + assert os.path.basename(subject_data[key]) == value + elif isinstance(value, list): + assert subject_data[key] is not None, f'Key {key} is None.' + assert len(subject_data[key]) == len(value) + for item, expected_item in zip(subject_data[key], value): + assert os.path.basename(item) == expected_item + else: + assert subject_data[key] is value diff --git a/src/fmripost_aroma/tests/utils.py b/src/fmripost_aroma/tests/utils.py new file mode 100644 index 0000000..00daf06 --- /dev/null +++ b/src/fmripost_aroma/tests/utils.py @@ -0,0 +1,12 @@ +"""Utility functions for tests.""" + +from pathlib import Path + + +def get_test_data_path(): + """Return the path to test datasets, terminated with separator. + + Test-related data are kept in tests folder in "data". + Based on function by Yaroslav Halchenko used in Neurosynth Python package. + """ + return Path(__file__).resolve().parent.parent.parent.parent / 'tests' / 'data' diff --git a/src/fmripost_aroma/utils/bids.py b/src/fmripost_aroma/utils/bids.py index 844fb7e..b263324 100644 --- a/src/fmripost_aroma/utils/bids.py +++ b/src/fmripost_aroma/utils/bids.py @@ -12,9 +12,19 @@ from fmripost_aroma.data import load as load_data -def extract_entities(file_list): +def extract_entities(file_list: str | list[str]) -> dict: """Return a dictionary of common entities given a list of files. + Parameters + ---------- + file_list : str | list[str] + File path or list of file paths. + + Returns + ------- + entities : dict + Dictionary of entities. + Examples -------- >>> extract_entities("sub-01/anat/sub-01_T1w.nii.gz") @@ -46,57 +56,106 @@ def _unique(inlist): def collect_derivatives( - raw_dir: Path | None, - derivatives_dir: Path, - entities: dict, + raw_dataset: Path | BIDSLayout | None, + derivatives_dataset: Path | BIDSLayout | None, + entities: dict | None, fieldmap_id: str | None, spec: dict | None = None, patterns: list[str] | None = None, -): - """Gather existing derivatives and compose a cache.""" + allow_multiple: bool = False, +) -> dict: + """Gather existing derivatives and compose a cache. + + Parameters + ---------- + raw_dataset : Path | BIDSLayout | None + Path to the raw dataset or a BIDSLayout instance. + derivatives_dataset : Path | BIDSLayout + Path to the derivatives dataset or a BIDSLayout instance. + entities : dict + Dictionary of entities to use for filtering. + fieldmap_id : str | None + Fieldmap ID to use for filtering. + spec : dict | None + Specification dictionary. + patterns : list[str] | None + List of patterns to use for filtering. + allow_multiple : bool + Allow multiple files to be returned for a given query. + + Returns + ------- + derivs_cache : dict + Dictionary with keys corresponding to the derivatives and values + corresponding to the file paths. + """ + if not entities: + entities = {} + if spec is None or patterns is None: - _spec, _patterns = tuple( - json.loads(load_data.readable('io_spec.json').read_text()).values() - ) + _spec = json.loads(load_data.readable('io_spec.json').read_text()) if spec is None: - spec = _spec + spec = _spec['queries'] + if patterns is None: - patterns = _patterns + patterns = _spec['patterns'] + # Search for derivatives data derivs_cache = defaultdict(list, {}) + if derivatives_dataset is not None: + layout = derivatives_dataset + if isinstance(derivatives_dataset, Path): + derivatives_dataset = BIDSLayout( + derivatives_dataset, + config=['bids', 'derivatives'], + validate=False, + ) - layout = BIDSLayout(derivatives_dir, config=['bids', 'derivatives'], validate=False) - derivatives_dir = Path(derivatives_dir) + for k, q in spec['derivatives'].items(): + # Combine entities with query. Query values override file entities. + query = {**entities, **q} + item = layout.get(return_type='filename', **query) + if not item: + derivs_cache[k] = None + elif not allow_multiple and len(item) > 1: + raise ValueError(f'Multiple files found for {k}: {item}') + else: + derivs_cache[k] = item[0] if len(item) == 1 else item - # Search for preprocessed BOLD data - for k, q in spec['baseline']['derivatives'].items(): - query = {**q, **entities} - item = layout.get(return_type='filename', **query) - if not item: - continue - derivs_cache[k] = item[0] if len(item) == 1 else item + for k, q in spec['transforms'].items(): + # Combine entities with query. Query values override file entities. + # TODO: Drop functional entities (task, run, etc.) from anat transforms. + query = {**entities, **q} + if k == 'boldref2fmap': + query['to'] = fieldmap_id + + item = layout.get(return_type='filename', **query) + if not item: + derivs_cache[k] = None + elif not allow_multiple and len(item) > 1: + raise ValueError(f'Multiple files found for {k}: {item}') + else: + derivs_cache[k] = item[0] if len(item) == 1 else item # Search for raw BOLD data - if not derivs_cache and raw_dir is not None: - raw_layout = BIDSLayout(raw_dir, config=['bids'], validate=False) - raw_dir = Path(raw_dir) + if not derivs_cache and raw_dataset is not None: + if isinstance(raw_dataset, Path): + raw_layout = BIDSLayout(raw_dataset, config=['bids'], validate=False) + else: + raw_layout = raw_dataset - for k, q in spec['baseline']['raw'].items(): - query = {**q, **entities} + for k, q in spec['raw'].items(): + # Combine entities with query. Query values override file entities. + query = {**entities, **q} item = raw_layout.get(return_type='filename', **query) if not item: - continue - derivs_cache[k] = item[0] if len(item) == 1 else item - - for xfm, q in spec['transforms'].items(): - query = {**q, **entities} - if xfm == 'boldref2fmap': - query['to'] = fieldmap_id - item = layout.get(return_type='filename', **q) - if not item: - continue - derivs_cache[xfm] = item[0] if len(item) == 1 else item + derivs_cache[k] = None + elif not allow_multiple and len(item) > 1: + raise ValueError(f'Multiple files found for {k}: {item}') + else: + derivs_cache[k] = item[0] if len(item) == 1 else item + return derivs_cache @@ -139,3 +198,172 @@ def collect_run_data( run_data[k] = layout.get_nearest(bold_file, **v) return run_data + + +def write_bidsignore(deriv_dir): + bids_ignore = ( + '*.html', + 'logs/', + 'figures/', # Reports + '*_xfm.*', # Unspecified transform files + '*.surf.gii', # Unspecified structural outputs + # Unspecified functional outputs + '*_boldref.nii.gz', + '*_bold.func.gii', + '*_mixing.tsv', + '*_timeseries.tsv', + ) + ignore_file = Path(deriv_dir) / '.bidsignore' + + ignore_file.write_text('\n'.join(bids_ignore) + '\n') + + +def write_derivative_description(bids_dir, deriv_dir, dataset_links=None): + import os + + from fmripost_aroma import __version__ + + DOWNLOAD_URL = f'https://github.com/nipreps/fmripost_aroma/archive/{__version__}.tar.gz' + + bids_dir = Path(bids_dir) + deriv_dir = Path(deriv_dir) + desc = { + 'Name': 'fMRIPost-AROMA- ICA-AROMA Postprocessing Outputs', + 'BIDSVersion': '1.9.0dev', + 'DatasetType': 'derivative', + 'GeneratedBy': [ + { + 'Name': 'fMRIPost-AROMA', + 'Version': __version__, + 'CodeURL': DOWNLOAD_URL, + } + ], + 'HowToAcknowledge': 'Please cite fMRIPost-AROMA when using these results.', + } + + # Keys that can only be set by environment + if 'FMRIPOST_AROMA_DOCKER_TAG' in os.environ: + desc['GeneratedBy'][0]['Container'] = { + 'Type': 'docker', + 'Tag': f"nipreps/fmriprep:{os.environ['FMRIPOST_AROMA__DOCKER_TAG']}", + } + if 'FMRIPOST_AROMA__SINGULARITY_URL' in os.environ: + desc['GeneratedBy'][0]['Container'] = { + 'Type': 'singularity', + 'URI': os.getenv('FMRIPOST_AROMA__SINGULARITY_URL'), + } + + # Keys deriving from source dataset + orig_desc = {} + fname = bids_dir / 'dataset_description.json' + if fname.exists(): + orig_desc = json.loads(fname.read_text()) + + if 'DatasetDOI' in orig_desc: + desc['SourceDatasets'] = [ + {'URL': f'https://doi.org/{orig_desc["DatasetDOI"]}', 'DOI': orig_desc['DatasetDOI']} + ] + if 'License' in orig_desc: + desc['License'] = orig_desc['License'] + + # Add DatasetLinks + if dataset_links: + desc['DatasetLinks'] = {k: str(v) for k, v in dataset_links.items()} + if 'templateflow' in dataset_links: + desc['DatasetLinks']['templateflow'] = 'https://github.com/templateflow/templateflow' + + Path.write_text(deriv_dir / 'dataset_description.json', json.dumps(desc, indent=4)) + + +def validate_input_dir(exec_env, bids_dir, participant_label, need_T1w=True): + # Ignore issues and warnings that should not influence FMRIPREP + import subprocess + import sys + import tempfile + + validator_config_dict = { + 'ignore': [ + 'EVENTS_COLUMN_ONSET', + 'EVENTS_COLUMN_DURATION', + 'TSV_EQUAL_ROWS', + 'TSV_EMPTY_CELL', + 'TSV_IMPROPER_NA', + 'VOLUME_COUNT_MISMATCH', + 'BVAL_MULTIPLE_ROWS', + 'BVEC_NUMBER_ROWS', + 'DWI_MISSING_BVAL', + 'INCONSISTENT_SUBJECTS', + 'INCONSISTENT_PARAMETERS', + 'BVEC_ROW_LENGTH', + 'B_FILE', + 'PARTICIPANT_ID_COLUMN', + 'PARTICIPANT_ID_MISMATCH', + 'TASK_NAME_MUST_DEFINE', + 'PHENOTYPE_SUBJECTS_MISSING', + 'STIMULUS_FILE_MISSING', + 'DWI_MISSING_BVEC', + 'EVENTS_TSV_MISSING', + 'TSV_IMPROPER_NA', + 'ACQTIME_FMT', + 'Participants age 89 or higher', + 'DATASET_DESCRIPTION_JSON_MISSING', + 'FILENAME_COLUMN', + 'WRONG_NEW_LINE', + 'MISSING_TSV_COLUMN_CHANNELS', + 'MISSING_TSV_COLUMN_IEEG_CHANNELS', + 'MISSING_TSV_COLUMN_IEEG_ELECTRODES', + 'UNUSED_STIMULUS', + 'CHANNELS_COLUMN_SFREQ', + 'CHANNELS_COLUMN_LOWCUT', + 'CHANNELS_COLUMN_HIGHCUT', + 'CHANNELS_COLUMN_NOTCH', + 'CUSTOM_COLUMN_WITHOUT_DESCRIPTION', + 'ACQTIME_FMT', + 'SUSPICIOUSLY_LONG_EVENT_DESIGN', + 'SUSPICIOUSLY_SHORT_EVENT_DESIGN', + 'MALFORMED_BVEC', + 'MALFORMED_BVAL', + 'MISSING_TSV_COLUMN_EEG_ELECTRODES', + 'MISSING_SESSION', + ], + 'error': ['NO_T1W'] if need_T1w else [], + 'ignoredFiles': ['/dataset_description.json', '/participants.tsv'], + } + # Limit validation only to data from requested participants + if participant_label: + all_subs = {s.name[4:] for s in bids_dir.glob('sub-*')} + selected_subs = {s[4:] if s.startswith('sub-') else s for s in participant_label} + bad_labels = selected_subs.difference(all_subs) + if bad_labels: + error_msg = ( + 'Data for requested participant(s) label(s) not found. Could ' + 'not find data for participant(s): %s. Please verify the requested ' + 'participant labels.' + ) + if exec_env == 'docker': + error_msg += ( + ' This error can be caused by the input data not being ' + 'accessible inside the docker container. Please make sure all ' + 'volumes are mounted properly (see https://docs.docker.com/' + 'engine/reference/commandline/run/#mount-volume--v---read-only)' + ) + if exec_env == 'singularity': + error_msg += ( + ' This error can be caused by the input data not being ' + 'accessible inside the singularity container. Please make sure ' + 'all paths are mapped properly (see https://www.sylabs.io/' + 'guides/3.0/user-guide/bind_paths_and_mounts.html)' + ) + raise RuntimeError(error_msg % ','.join(bad_labels)) + + ignored_subs = all_subs.difference(selected_subs) + if ignored_subs: + for sub in ignored_subs: + validator_config_dict['ignoredFiles'].append(f'/sub-{sub}/**') + with tempfile.NamedTemporaryFile(mode='w+', suffix='.json') as temp: + temp.write(json.dumps(validator_config_dict)) + temp.flush() + try: + subprocess.check_call(['bids-validator', str(bids_dir), '-c', temp.name]) # noqa: S607 + except FileNotFoundError: + print('bids-validator does not appear to be installed', file=sys.stderr) diff --git a/src/fmripost_aroma/utils/utils.py b/src/fmripost_aroma/utils/utils.py index 4836a59..dcaf7c3 100644 --- a/src/fmripost_aroma/utils/utils.py +++ b/src/fmripost_aroma/utils/utils.py @@ -441,3 +441,29 @@ def _get_wf_name(bold_fname, prefix): fname = split_filename(bold_fname)[1] fname_nosub = '_'.join(fname.split('_')[1:-1]) return f"{prefix}_{fname_nosub.replace('-', '_')}_wf" + + +def update_dict(orig_dict, new_dict): + """Update dictionary with values from another dictionary. + + Parameters + ---------- + orig_dict : dict + Original dictionary. + new_dict : dict + Dictionary with new values. + + Returns + ------- + updated_dict : dict + Updated dictionary. + """ + updated_dict = orig_dict.copy() + for key, value in new_dict.items(): + if (orig_dict.get(key) is not None) and (value is not None): + print(f'Updating {key} from {orig_dict[key]} to {value}') + updated_dict[key].update(value) + elif value is not None: + updated_dict[key] = value + + return updated_dict diff --git a/src/fmripost_aroma/workflows/aroma.py b/src/fmripost_aroma/workflows/aroma.py index 8c623bc..b048f2c 100644 --- a/src/fmripost_aroma/workflows/aroma.py +++ b/src/fmripost_aroma/workflows/aroma.py @@ -114,7 +114,7 @@ def init_ica_aroma_wf( workflow = Workflow(name=_get_wf_name(bold_file, 'aroma')) workflow.__postdesc__ = f"""\ Automatic removal of motion artifacts using independent component analysis -[ICA-AROMA, @aroma] was performed on the *preprocessed BOLD on MNI152NLin6Asym space* +[ICA-AROMA, @ica_aroma] was performed on the *preprocessed BOLD on MNI152NLin6Asym space* time-series after removal of non-steady state volumes and spatial smoothing with a nonlinear filter that preserves underlying structure [SUSAN, @susan], using a FWHM of {susan_fwhm} mm. @@ -259,6 +259,7 @@ def init_ica_aroma_wf( ds_report_ica_aroma = pe.Node( DerivativesDataSink( + base_directory=config.execution.fmripost_aroma_dir, source_file=bold_file, datatype='figures', desc='aroma', @@ -289,6 +290,7 @@ def init_ica_aroma_wf( ds_components = pe.Node( DerivativesDataSink( + base_directory=config.execution.fmripost_aroma_dir, source_file=bold_file, datatype='func', desc='melodic', @@ -304,6 +306,7 @@ def init_ica_aroma_wf( ds_mixing = pe.Node( DerivativesDataSink( + base_directory=config.execution.fmripost_aroma_dir, source_file=bold_file, datatype='func', desc='melodic', @@ -319,6 +322,7 @@ def init_ica_aroma_wf( ds_aroma_features = pe.Node( DerivativesDataSink( + base_directory=config.execution.fmripost_aroma_dir, source_file=bold_file, datatype='func', desc='aroma', @@ -339,6 +343,7 @@ def init_ica_aroma_wf( ds_aroma_confounds = pe.Node( DerivativesDataSink( + base_directory=config.execution.fmripost_aroma_dir, source_file=bold_file, datatype='func', desc='melodic', @@ -374,7 +379,8 @@ def init_denoise_wf(bold_file): fields=[ 'bold_file', 'bold_mask', - 'confounds', + 'mixing', + 'classifications', 'skip_vols', 'spatial_reference', ], @@ -400,8 +406,10 @@ def init_denoise_wf(bold_file): ) workflow.connect([ (inputnode, denoise, [ - ('confounds', 'confounds'), + ('mixing', 'mixing'), + ('classifications', 'metrics'), ('bold_mask', 'mask_file'), + ('skip_vols', 'skip_vols'), ]), (rm_non_steady_state, denoise, [('bold_cut', 'bold_file')]), ]) # fmt:skip @@ -420,6 +428,7 @@ def init_denoise_wf(bold_file): ds_denoised = pe.Node( DerivativesDataSink( + base_directory=config.execution.fmripost_aroma_dir, source_file=bold_file, desc=f'{denoise_method}Denoised', ), @@ -430,7 +439,7 @@ def init_denoise_wf(bold_file): workflow.connect([ # spatial_reference needs to be parsed into space, cohort, res, den, etc. (inputnode, ds_denoised, [('spatial_reference', 'space')]), - (add_non_steady_state, ds_denoised, [('bold_add', 'denoised_file')]), + (add_non_steady_state, ds_denoised, [('bold_add', 'in_file')]), ]) # fmt:skip return workflow diff --git a/src/fmripost_aroma/workflows/base.py b/src/fmripost_aroma/workflows/base.py index d5eb1da..5168b28 100644 --- a/src/fmripost_aroma/workflows/base.py +++ b/src/fmripost_aroma/workflows/base.py @@ -29,18 +29,15 @@ """ -import sys -import warnings +import os from copy import deepcopy -from nipype.interfaces import utility as niu +import yaml from nipype.pipeline import engine as pe from packaging.version import Version from fmripost_aroma import config -from fmripost_aroma.interfaces.bids import DerivativesDataSink -from fmripost_aroma.interfaces.reportlets import AboutSummary, SubjectSummary -from fmripost_aroma.utils.utils import _get_wf_name +from fmripost_aroma.utils.utils import _get_wf_name, update_dict from fmripost_aroma.workflows.resampling import init_resample_volumetric_wf @@ -118,11 +115,6 @@ def init_single_subject_wf(subject_id: str): subject_id : :obj:`str` Subject label for this single-subject workflow. - Inputs - ------ - subjects_dir : :obj:`str` - FreeSurfer's ``$SUBJECTS_DIR``. - Notes ----- 1. Load fMRIPost-AROMA config file. @@ -146,17 +138,14 @@ def init_single_subject_wf(subject_id: str): 8. Warp BOLD to requested output spaces and denoise with ICA-AROMA. """ + from bids.utils import listify from niworkflows.engine.workflows import LiterateWorkflow as Workflow - from niworkflows.interfaces.bids import BIDSDataGrabber, BIDSInfo from niworkflows.interfaces.nilearn import NILEARN_VERSION - from niworkflows.utils.misc import fix_multi_T1w_source_name from niworkflows.utils.spaces import Reference from fmripost_aroma.utils.bids import collect_derivatives, extract_entities from fmripost_aroma.workflows.aroma import init_denoise_wf, init_ica_aroma_wf - spaces = config.workflow.spaces - workflow = Workflow(name=f'sub_{subject_id}_wf') workflow.__desc__ = f""" Results included in this manuscript come from postprocessing @@ -187,10 +176,29 @@ def init_single_subject_wf(subject_id: str): """ - subject_data = collect_derivatives( - raw_dir=config.execution.layout, - entities=config.execution.bids_filters, - ) + if config.execution.derivatives: + # Raw dataset + derivatives dataset + config.loggers.workflow.info('Raw+derivatives workflow mode enabled') + subject_data = collect_derivatives( + raw_dataset=config.execution.layout, + derivatives_dataset=None, + entities=config.execution.bids_filters, + fieldmap_id=None, + allow_multiple=True, + ) + subject_data['bold'] = listify(subject_data['bold_raw']) + else: + # Derivatives dataset only + config.loggers.workflow.info('Derivatives-only workflow mode enabled') + subject_data = collect_derivatives( + raw_dataset=None, + derivatives_dataset=config.execution.layout, + entities=config.execution.bids_filters, + fieldmap_id=None, + allow_multiple=True, + ) + # Patch standard-space BOLD files into 'bold' key + subject_data['bold'] = listify(subject_data['bold_std']) # Make sure we always go through these two checks if not subject_data['bold']: @@ -198,82 +206,10 @@ def init_single_subject_wf(subject_id: str): raise RuntimeError( f"No BOLD images found for participant {subject_id} and " f"task {task_id if task_id else ''}. " - "All workflows require BOLD images." + "All workflows require BOLD images. " + f"Please check your BIDS filters: {config.execution.bids_filters}." ) - if subject_data['roi']: - warnings.warn( - f"Lesion mask {subject_data['roi']} found. " - "Future versions of fMRIPost-AROMA will use alternative conventions. " - "Please refer to the documentation before upgrading.", - FutureWarning, - stacklevel=1, - ) - - inputnode = pe.Node(niu.IdentityInterface(fields=['subjects_dir']), name='inputnode') - - bidssrc = pe.Node( - BIDSDataGrabber( - subject_data=subject_data, - subject_id=subject_id, - ), - name='bidssrc', - ) - - bids_info = pe.Node( - BIDSInfo(bids_dir=config.execution.bids_dir, bids_validate=False), - name='bids_info', - ) - - summary = pe.Node( - SubjectSummary( - std_spaces=['MNI152NLin6Asym'], - nstd_spaces=None, - ), - name='summary', - run_without_submitting=True, - ) - - about = pe.Node( - AboutSummary(version=config.environment.version, command=' '.join(sys.argv)), - name='about', - run_without_submitting=True, - ) - - ds_report_summary = pe.Node( - DerivativesDataSink( - base_directory=config.execution.fmripost_aroma_dir, - desc='summary', - datatype='figures', - dismiss_entities=('echo',), - ), - name='ds_report_summary', - run_without_submitting=True, - ) - - ds_report_about = pe.Node( - DerivativesDataSink( - base_directory=config.execution.fmripost_aroma_dir, - desc='about', - datatype='figures', - dismiss_entities=('echo',), - ), - name='ds_report_about', - run_without_submitting=True, - ) - - workflow.connect([ - (bidssrc, bids_info, [(('t1w', fix_multi_T1w_source_name), 'in_file')]), - # Reporting connections - (inputnode, summary, [('subjects_dir', 'subjects_dir')]), - (bidssrc, summary, [('t1w', 't1w'), ('t2w', 't2w'), ('bold', 'bold')]), - (bids_info, summary, [('subject', 'subject_id')]), - (bidssrc, ds_report_summary, [(('t1w', fix_multi_T1w_source_name), 'source_file')]), - (bidssrc, ds_report_about, [(('t1w', fix_multi_T1w_source_name), 'source_file')]), - (summary, ds_report_summary, [('out_report', 'in_file')]), - (about, ds_report_about, [('out_report', 'in_file')]), - ]) # fmt:skip - # Append the functional section to the existing anatomical excerpt # That way we do not need to stream down the number of bold datasets func_pre_desc = f""" @@ -284,7 +220,8 @@ def init_single_subject_wf(subject_id: str): """ for bold_file in subject_data['bold']: - ica_aroma_wf = init_ica_aroma_wf(bold_file=bold_file) + bold_metadata = config.execution.layout.get_metadata(bold_file) + ica_aroma_wf = init_ica_aroma_wf(bold_file=bold_file, metadata=bold_metadata) ica_aroma_wf.__desc__ = func_pre_desc + (ica_aroma_wf.__desc__ or '') entities = extract_entities(bold_file) @@ -292,14 +229,29 @@ def init_single_subject_wf(subject_id: str): functional_cache = {} if config.execution.derivatives: # Collect native-space derivatives and transforms - for deriv_dir in config.execution.derivatives: - functional_cache.update( + for deriv_dir in config.execution.derivatives.values(): + functional_cache = update_dict( + functional_cache, collect_derivatives( - derivatives_dir=deriv_dir, + raw_dataset=None, + derivatives_dataset=deriv_dir, entities=entities, + fieldmap_id=None, + allow_multiple=False, ), ) + if not functional_cache['confounds']: + if config.workflow.dummy_scans is None: + raise ValueError( + 'No confounds detected. ' + 'Automatical dummy scan detection cannot be performed. ' + 'Please set the `--dummy-scans` flag explicitly.' + ) + + # TODO: Calculate motion parameters from motion correction transform + raise ValueError('Motion parameters cannot be extracted from transforms yet.') + # Resample to MNI152NLin6Asym:res-2, for ICA-AROMA classification resample_raw_wf = init_resample_volumetric_wf( bold_file=bold_file, @@ -319,44 +271,65 @@ def init_single_subject_wf(subject_id: str): # Only derivatives dataset was passed in, so we expected standard-space derivatives functional_cache.update( collect_derivatives( - derivatives_dir=config.execution.layout, + raw_dataset=None, + derivatives_dataset=config.execution.layout, entities=entities, + fieldmap_id=None, + allow_multiple=False, ), ) ica_aroma_wf.inputs.inputnode.bold_std = functional_cache['bold_std'] ica_aroma_wf.inputs.inputnode.bold_mask_std = functional_cache['bold_mask_std'] workflow.add_nodes([ica_aroma_wf]) - functional_cache['skip_vols'] = ( - config.workflow.dummy_scans or functional_cache['skip_vols'] + config.loggers.workflow.info( + ( + f'Collected run data for {os.path.basename(bold_file)}:\n' + f'{yaml.dump(functional_cache, default_flow_style=False, indent=4)}' + ), ) - ica_aroma_wf.inputs.inputnode.confounds = functional_cache['confounds'] - ica_aroma_wf.inputs.inputnode.skip_vols = functional_cache['skip_vols'] - if config.workflow.denoise_method: - for space in spaces: - # Warp each BOLD run to requested output spaces - resample_to_space_wf = init_resample_volumetric_wf( - bold_file=bold_file, - functional_cache=functional_cache, - space=space, - name=_get_wf_name(bold_file, f'resample_{space}'), + if config.workflow.dummy_scans is not None: + skip_vols = config.workflow.dummy_scans + else: + if not functional_cache['confounds']: + raise ValueError( + 'No confounds detected. ' + 'Automatical dummy scan detection cannot be performed. ' + 'Please set the `--dummy-scans` flag explicitly.' ) + skip_vols = get_nss(functional_cache['confounds']) - # Now denoise the output-space BOLD data using ICA-AROMA - denoise_wf = init_denoise_wf(bold_file=bold_file) - denoise_wf.inputs.inputnode.bold_mask = functional_cache['bold_mask'] - denoise_wf.inputs.inputnode.skip_vols = functional_cache['skip_vols'] - workflow.connect([ - (resample_to_space_wf, denoise_wf, [ - ('bold_std', 'inputnode.bold_file'), - ('bold_mask_std', 'inputnode.bold_mask'), - ('spatial_reference', 'inputnode.spatial_reference'), - ]), - (ica_aroma_wf, denoise_wf, [ - ('outputnode.aroma_confounds', 'inputnode.confounds'), - ]), - ]) # fmt:skip + ica_aroma_wf.inputs.inputnode.confounds = functional_cache['confounds'] + ica_aroma_wf.inputs.inputnode.skip_vols = skip_vols + + if config.workflow.denoise_method: + # for space in spaces: + # # Warp each BOLD run to requested output spaces + # resample_to_space_wf = init_resample_volumetric_wf( + # bold_file=bold_file, + # functional_cache=functional_cache, + # space=space, + # name=_get_wf_name(bold_file, f'resample_{space}'), + # ) + + # Now denoise the output-space BOLD data using ICA-AROMA + denoise_wf = init_denoise_wf(bold_file=bold_file) + denoise_wf.inputs.inputnode.skip_vols = skip_vols + denoise_wf.inputs.inputnode.bold_file = functional_cache['bold_std'] + denoise_wf.inputs.inputnode.bold_mask = functional_cache['bold_mask_std'] + denoise_wf.inputs.inputnode.spatial_reference = 'MNI152NLin6Asym' + workflow.connect([ + # (resample_to_space_wf, denoise_wf, [ + # ('bold_std', 'inputnode.bold_file'), + # ('bold_mask_std', 'inputnode.bold_mask'), + # ('spatial_reference', 'inputnode.spatial_reference'), + # ]), + (ica_aroma_wf, denoise_wf, [ + ('outputnode.mixing', 'inputnode.mixing'), + ('outputnode.aroma_features', 'inputnode.classifications'), + ]), + ]) # fmt:skip return clean_datasinks(workflow) @@ -371,3 +344,24 @@ def clean_datasinks(workflow: pe.Workflow) -> pe.Workflow: if node.split('.')[-1].startswith('ds_'): workflow.get_node(node).interface.out_path_base = '' return workflow + + +def get_nss(confounds_file): + """Get number of non-steady state volumes.""" + import numpy as np + import pandas as pd + + df = pd.read_table(confounds_file) + + nss_cols = [c for c in df.columns if c.startswith('non_steady_state_outlier')] + + dummy_scans = 0 + if nss_cols: + initial_volumes_df = df[nss_cols] + dummy_scans = np.any(initial_volumes_df.to_numpy(), axis=1) + dummy_scans = np.where(dummy_scans)[0] + + # reasonably assumes all NSS volumes are contiguous + dummy_scans = int(dummy_scans[-1] + 1) + + return dummy_scans