Skip to content
New issue

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

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

Already on GitHub? Sign in to your account

[ENH] integrate load_confounds into first_level_from_bids #4103

Merged
merged 21 commits into from
Dec 18, 2023
Merged
Show file tree
Hide file tree
Changes from 4 commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
2 changes: 1 addition & 1 deletion .flake8
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,7 @@ docstring-convention = numpy
max-line-length = 79
max_complexity = 43
max_function_length = 407
max_parameters_amount = 26
max_parameters_amount = 27
max_returns_amount = 10
# For PEP8 error codes see
# http://pep8.readthedocs.org/en/latest/intro.html#error-codes
Expand Down
1 change: 0 additions & 1 deletion .pre-commit-config.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,6 @@ exclude: |
| nilearn/interfaces/fmriprep/load_confounds_compcor.py
| nilearn/interfaces/fmriprep/tests/test_load_confounds_utils.py
| nilearn/interfaces/fmriprep/tests/test_load_confounds_strategy.py
| nilearn/interfaces/fmriprep/tests/utils.py
| nilearn/interfaces/fmriprep/tests/test_load_confounds.py
| nilearn/maskers/nifti_spheres_masker.py
| nilearn/maskers/nifti_labels_masker.py
Expand Down
39 changes: 39 additions & 0 deletions nilearn/_utils/bids.py
Copy link
Collaborator Author

Choose a reason for hiding this comment

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

had to move the create_bids_filename in a different module to help avoid circular imports

Original file line number Diff line number Diff line change
@@ -0,0 +1,39 @@
from nilearn.interfaces.bids._utils import _bids_entities


def create_bids_filename(fields, entities_to_include=None):
"""Create BIDS filename from dictionary of entity-label pairs.

Parameters
----------
fields : :obj:`dict` of :obj:`str`
Dictionary of entity-label pairs, for example:

{
"suffix": "T1w",
"extension": "nii.gz",
"entities": {"acq": "ap",
"desc": "preproc"}
}.

Returns
-------
BIDS filename : :obj:`str`

"""
if entities_to_include is None:
entities_to_include = _bids_entities()["raw"]

filename = ""

for key in entities_to_include:
if key in fields["entities"]:
value = fields["entities"][key]
if value not in (None, ""):
filename += f"{key}-{value}_"
if "suffix" in fields:
filename += f"{fields['suffix']}"
if "extension" in fields:
filename += f".{fields['extension']}"

return filename
77 changes: 24 additions & 53 deletions nilearn/_utils/data_gen.py
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,9 @@

from nilearn import datasets, image, maskers, masking
from nilearn._utils import as_ndarray, logger
from nilearn._utils.bids import create_bids_filename
from nilearn.interfaces.bids._utils import _bids_entities, _check_bids_label
from nilearn.interfaces.fmriprep.tests._testing import get_legal_confound


def generate_mni_space_img(n_scans=1, res=30, random_state=0, mask_dilation=2):
Expand Down Expand Up @@ -735,9 +737,15 @@ def _basic_confounds(length, random_state=0):

"""
rand_gen = check_random_state(random_state)
columns = ['csf', 'white_matter', 'global_signal',
'rot_x', 'rot_y', 'rot_z',
'trans_x', 'trans_y', 'trans_z']
columns = ['csf',
'white_matter',
'global_signal',
'rot_x',
'rot_y',
'rot_z',
'trans_x',
'trans_y',
'trans_z']
data = rand_gen.rand(length, len(columns))
confounds = pd.DataFrame(data, columns=columns)
return confounds
Expand Down Expand Up @@ -1198,46 +1206,6 @@ def _listify(n):
return [""] if n <= 0 else [f"{label:02}" for label in range(1, n + 1)]


def _create_bids_filename(
fields, entities_to_include=None
):
"""Create BIDS filename from dictionary of entity-label pairs.

Parameters
----------
fields : :obj:`dict` of :obj:`str`
Dictionary of entity-label pairs, for example:

{
"suffix": "T1w",
"extension": "nii.gz",
"entities": {"acq": "ap",
"desc": "preproc"}
}.

Returns
-------
BIDS filename : :obj:`str`

"""
if entities_to_include is None:
entities_to_include = _bids_entities()["raw"]

filename = ""

for key in entities_to_include:
if key in fields["entities"]:
value = fields["entities"][key]
if value not in (None, ""):
filename += f"{key}-{value}_"
if "suffix" in fields:
filename += f"{fields['suffix']}"
if "extension" in fields:
filename += f".{fields['extension']}"

return filename


def _init_fields(subject,
session,
task,
Expand Down Expand Up @@ -1265,7 +1233,7 @@ def _init_fields(subject,

See Also
--------
_create_bids_filename
create_bids_filename

"""
fields = {
Expand Down Expand Up @@ -1302,7 +1270,7 @@ def _write_bids_raw_anat(subses_dir, subject, session) -> None:
"extension": "nii.gz",
"entities": {"sub": subject, "ses": session},
}
(anat_path / _create_bids_filename(fields)).write_text("")
(anat_path / create_bids_filename(fields)).write_text("")


def _write_bids_raw_func(
Expand Down Expand Up @@ -1330,7 +1298,7 @@ def _write_bids_raw_func(

"""
n_time_points = 100
bold_path = func_path / _create_bids_filename(fields)
bold_path = func_path / create_bids_filename(fields)
write_fake_bold_img(
bold_path,
[n_voxels, n_voxels, n_voxels, n_time_points],
Expand All @@ -1339,12 +1307,12 @@ def _write_bids_raw_func(

repetition_time = 1.5
fields["extension"] = "json"
param_path = func_path / _create_bids_filename(fields)
param_path = func_path / create_bids_filename(fields)
param_path.write_text(json.dumps({"RepetitionTime": repetition_time}))

fields["suffix"] = "events"
fields["extension"] = "tsv"
events_path = func_path / _create_bids_filename(fields)
events_path = func_path / create_bids_filename(fields)
basic_paradigm().to_csv(events_path, sep="\t", index=None)


Expand Down Expand Up @@ -1390,12 +1358,15 @@ def _write_bids_derivative_func(
if confounds_tag is not None:
fields["suffix"] = confounds_tag
fields["extension"] = "tsv"
confounds_path = func_path / _create_bids_filename(
confounds_path = func_path / create_bids_filename(
fields=fields, entities_to_include=_bids_entities()["raw"]
)
_basic_confounds(length=n_time_points, random_state=rand_gen).to_csv(
confounds_path, sep="\t", index=None
confounds, metadata = get_legal_confound()
confounds.to_csv(
confounds_path, sep="\t", index=None, encoding="utf-8"
)
with open(confounds_path.with_suffix(".json"), "w") as f:
json.dump(metadata, f)
Comment on lines +1367 to +1372
Copy link
Collaborator Author

Choose a reason for hiding this comment

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

One issue with this approach is that the "legal_confounds" have a set number of time points so when creating a fake bids dataset, we end up with images that have a number of time points that does not match the number of time points in the confounds.

This does not affect any tests AFIACT but this may lead to confusing errors when testing down the line.

Copy link
Member

Choose a reason for hiding this comment

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

You mean, because of scrubbing ? Sorry if I miss something obvious.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

nope

let me try to rephrase

the way to generate "fake confounds" for the fake fmriprep datasets we use for testing would only create 6 confounds for the realignment parameters filled with random data and for a specified number of time points

to allow testing the load_confounds we need more realistic confounds with more columns with names that match what's in an actual fmriprep dataset

to do this I reuse the strategy used to test the load_confounds functions: use an actual confound file from an fmriprep dataset and copy its content every time it is needed in the fake fmriprep dataset

but this "template" confound file has only a limited number of time points

so we end up with fake fmriprep datasets that have nifti images with 100 volumes but with confounds with only 30 time points

possible solutions:

  • easy: set the number of volumes to match the number of time points in the confounds
  • hard(er): adapt the content of the confounds to the number of time points

hope this is clearer

for now I will go for the easy solution but we may have to implement the harder option in the future if we want to test more "exotic" stuff

Copy link
Member

Choose a reason for hiding this comment

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

Let's go for the easy one. The number of volumes should be a parameters of the data simulation function anyhow ?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

for now it is not: it is hard coded. I would keep it that way for now until we need more flexibility during testing.

but I will change the place where it is hard coded so it is easier to adapt in the future. will also add a comment to explain why this value was chosen.

Comment on lines +1371 to +1372
Copy link
Collaborator Author

Choose a reason for hiding this comment

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

minor change: using legal_confounds allows to add metada files for the confounds in the fake bids derivatives. some tests had to be changed to account for this.


fields["suffix"] = "bold"
fields["extension"] = "nii.gz"
Expand All @@ -1416,7 +1387,7 @@ def _write_bids_derivative_func(
fields["entities"]["space"] = space
fields["entities"]["desc"] = desc

bold_path = func_path / _create_bids_filename(
bold_path = func_path / create_bids_filename(
fields=fields, entities_to_include=entities_to_include
)
write_fake_bold_img(bold_path, shape=shape, random_state=rand_gen)
Expand All @@ -1426,7 +1397,7 @@ def _write_bids_derivative_func(
fields["entities"].pop("desc")
for hemi in ["L", "R"]:
fields["entities"]["hemi"] = hemi
gifti_path = func_path / _create_bids_filename(
gifti_path = func_path / create_bids_filename(
fields=fields,
entities_to_include=entities_to_include
)
Expand Down
10 changes: 5 additions & 5 deletions nilearn/_utils/tests/test_data_gen.py
Original file line number Diff line number Diff line change
Expand Up @@ -211,8 +211,8 @@ def test_fake_bids_derivatives_with_session_and_runs(
)

all_files = list(bids_path.glob("derivatives/sub-*/ses-*/*/*"))
# per subject: (1 confound + 3 bold + 2 gifti) per run per session
n_derivatives_files_expected = n_sub * (6 * sum(n_runs) * n_ses)
# per subject: (2 confound + 3 bold + 2 gifti) per run per session
n_derivatives_files_expected = n_sub * (7 * sum(n_runs) * n_ses)
assert len(all_files) == n_derivatives_files_expected


Expand Down Expand Up @@ -375,10 +375,10 @@ def test_fake_bids_extra_raw_entity(tmp_path):
)

all_files = list(bids_path.glob("derivatives/sub-*/ses-*/*/*"))
# per subject: (1 confound + 3 bold + 2 gifti)
# per subject: (2 confound + 3 bold + 2 gifti)
# per run per session per entity
n_derivatives_files_expected = (
n_sub * (6 * sum(n_runs) * n_ses) * len(entities["acq"])
n_sub * (7 * sum(n_runs) * n_ses) * len(entities["acq"])
)
assert len(all_files) == n_derivatives_files_expected

Expand Down Expand Up @@ -420,7 +420,7 @@ def test_fake_bids_extra_derivative_entity(tmp_path):
# 1 confound per run per session
# + (3 bold + 2 gifti) per run per session per entity
n_derivatives_files_expected = n_sub * (
1 * sum(n_runs) * n_ses
2 * sum(n_runs) * n_ses
+ 5 * sum(n_runs) * n_ses * len(entities["res"])
)
assert len(all_files) == n_derivatives_files_expected
Expand Down
58 changes: 51 additions & 7 deletions nilearn/glm/first_level/first_level.py
Original file line number Diff line number Diff line change
Expand Up @@ -48,6 +48,7 @@
_infer_repetition_time_from_dataset,
_infer_slice_timing_start_time_from_dataset,
)
from nilearn.interfaces.fmriprep.load_confounds import load_confounds


def mean_scaling(Y, axis=0):
Expand Down Expand Up @@ -957,6 +958,7 @@
n_jobs=1,
minimize_memory=True,
derivatives_folder="derivatives",
**kwargs,
):
"""Create FirstLevelModel objects and fit arguments \
from a :term:`BIDS` dataset.
Expand Down Expand Up @@ -1050,6 +1052,11 @@
derivatives_folder=derivatives_folder,
)

kwargs_load_confounds = _check_kwargs_load_confounds(**kwargs)

# TODO check compatibility of kwargs with FirstLevelModel
# for example high_pass and strategy high_pass

derivatives_path = Path(dataset_path) / derivatives_folder

# Get metadata for models.
Expand Down Expand Up @@ -1216,11 +1223,8 @@
img_filters=img_filters,
imgs=imgs,
verbose=verbose,
kwargs_load_confounds=kwargs_load_confounds,
)
if confounds:
confounds = [
pd.read_csv(c, sep="\t", index_col=None) for c in confounds
]
models_confounds.append(confounds)

return models, models_run_imgs, models_events, models_confounds
Expand Down Expand Up @@ -1433,6 +1437,7 @@
img_filters,
imgs,
verbose,
kwargs_load_confounds,
):
"""Get confounds.tsv files for a given subject, task and filters.

Expand Down Expand Up @@ -1462,8 +1467,7 @@

Returns
-------
confounds : :obj:`list` of :obj:`str` or None
List of fullpath to the confounds.tsv files
confounds : :obj:`list` of :class:`pandas.DataFrame`

"""
filters = _make_bids_files_filter(
Expand All @@ -1488,7 +1492,17 @@
filters=filters,
)
_check_confounds_list(confounds=confounds, imgs=imgs)
return confounds or None

if confounds:
if kwargs_load_confounds is None:
confounds = [
pd.read_csv(c, sep="\t", index_col=None) for c in confounds
]
return confounds or None

confounds, _ = load_confounds(img_files=imgs, **kwargs_load_confounds)

return confounds


def _check_confounds_list(confounds, imgs):
Expand Down Expand Up @@ -1613,6 +1627,36 @@
_check_bids_label(filter_[1])


def _check_kwargs_load_confounds(**kwargs):
# reuse the default from nilearn.interface.fmriprep.load_confounds
defaults = {
"motion": "full",
"scrub": 5,
"fd_threshold": 0.2,
"std_dvars_threshold": 3,
"wm_csf": "basic",
"global_signal": "basic",
"compcor": "anat_combined",
"n_compcor": "all",
"ica_aroma": "full",
"demean": True,
}

if "confounds_strategy" not in kwargs:
return None
if kwargs["confounds_strategy"] is None:
return None

Check warning on line 1648 in nilearn/glm/first_level/first_level.py

View check run for this annotation

Codecov / codecov/patch

nilearn/glm/first_level/first_level.py#L1648

Added line #L1648 was not covered by tests
Remi-Gau marked this conversation as resolved.
Show resolved Hide resolved

kwargs_load_confounds = {
key: defaults[key]
if f"confounds_{key}" not in kwargs
else kwargs[f"confounds_{key}"]
for key in defaults
}

return kwargs_load_confounds


def _make_bids_files_filter(
task_label,
space_label=None,
Expand Down
27 changes: 27 additions & 0 deletions nilearn/glm/tests/test_first_level.py
Original file line number Diff line number Diff line change
Expand Up @@ -1860,3 +1860,30 @@ def test_missing_trial_type_column_warning(tmp_path_factory):
"No column named 'trial_type' found" in r.message.args[0]
for r in record
)


def test_first_level_from_bids_load_confounds(tmp_path):
"""Test several BIDS structure."""
n_sub = 2

bids_path = create_fake_bids_dataset(
base_dir=tmp_path, n_sub=n_sub, n_ses=2, tasks=["main"], n_runs=[2]
)

models, m_imgs, m_events, m_confounds = first_level_from_bids(
dataset_path=bids_path,
task_label="main",
space_label="MNI",
img_filters=[("desc", "preproc")],
slice_time_ref=None,
confounds_strategy=("motion", "wm_csf", "scrub"),
confounds_motion="full",
confounds_wm_csf="basic",
confounds_scrub=1,
confounds_fd_threshold=0.2,
confounds_std_dvars_threshold=3,
)
Copy link
Collaborator Author

@Remi-Gau Remi-Gau Nov 9, 2023

Choose a reason for hiding this comment

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

this what the "API" could look like to select some confounds: does it look sensible?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

@htwangtw @ymzayek @bthirion
keeping this as a draft so you can discuss API, implementation...

Will update the doc once the dust settles.

Copy link
Member

Choose a reason for hiding this comment

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

I guess this could be made lighter if we expect to have predefined configurations for confounds, but I don't think that we're at that point.
At least I find the current API quite explicit.

Copy link
Member

Choose a reason for hiding this comment

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

Actually, we might get rid of confounds_strategy, because the next three arguments are redundant with the provided list ?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Note this test just show a subset of the possible arguments that can be passed to load_confounds

I guess this could be made lighter if we expect to have predefined configurations for confounds, but I don't think that we're at that point.

Actually, we might get rid of confounds_strategy, because the next three arguments are redundant with the provided list ?

Actually I am just reusing the API from load_confounds, so it can almost be passed as is and we can let load_confounds do the argument validation

https://nilearn.github.io/dev/modules/generated/nilearn.interfaces.fmriprep.load_confounds.html

In short strategy defines what type of confounds to include and all the other parameters give more details on how to include them.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

will update the doc string to try to explain this so we can see if it makes sense


_check_output_first_level_from_bids(
n_sub, models, m_imgs, m_events, m_confounds
)