Skip to content

[MRG] FIX: Output NaN columns for CompCor on failure #1443

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

Merged
merged 13 commits into from
Dec 18, 2018
Merged
2 changes: 1 addition & 1 deletion docs/environment.yml
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@ dependencies:
- python-dateutil
- pydot>=1.2.3
- cython
- nipype>=1.1.6
- nipype>=1.1.7

- pip:
- sphinx-argparse
Expand Down
2 changes: 1 addition & 1 deletion fmriprep/__about__.py
Original file line number Diff line number Diff line change
Expand Up @@ -88,7 +88,7 @@
'indexed_gzip>=0.8.8',
'nibabel>=2.2.1',
'nilearn',
'nipype>=1.1.6',
'nipype>=1.1.7',
'nitime',
'niworkflows>=0.5.1,<0.5.2',
'numpy',
Expand Down
111 changes: 83 additions & 28 deletions fmriprep/interfaces/patches.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,49 +12,104 @@

from numpy.linalg.linalg import LinAlgError
from nipype.algorithms import confounds as nac
from nipype.interfaces.base import File
from nipype.interfaces.mixins import reporting


class RobustACompCor(nac.ACompCor):
"""
Runs aCompCor several times if it suddenly fails with
https://github.com/poldracklab/fmriprep/issues/776
class RetryCompCorInputSpecMixin(reporting.ReportCapableInputSpec):
out_report = File('report.html', usedefault=True, hash_files=False,
desc='filename for warning HTML snippet')

"""

class RetryCompCorMixin(reporting.ReportCapableInterface):
def _run_interface(self, runtime):
warn = self.inputs.failure_mode == 'NaN'

failures = 0
save_exc = None
while True:
success = True
# Identifiy success/failure in both error and NaN mode
try:
runtime = super(RobustACompCor, self)._run_interface(runtime)
runtime = super()._run_interface(runtime)
if warn and self._is_allnans():
success = False
except LinAlgError as exc:
success = False
save_exc = exc

if success:
break
except LinAlgError:
failures += 1
if failures > 10:
raise
start = (failures - 1) * 10
sleep(randint(start + 4, start + 10))

failures += 1
if failures > 10:
if warn:
break
raise save_exc
start = (failures - 1) * 10
sleep(randint(start + 4, start + 10))

return runtime

def _is_allnans(self):
import numpy as np
outputs = self._list_outputs()
components = np.loadtxt(outputs['components_file'], skiprows=1)
return np.isnan(components).all()

def _generate_report(self):
snippet = '<!-- {} completed without error -->'.format(self._header)
if self._is_allnans():
snippet = '''\
<p class="elem-desc">
Warning: {} components could not be estimated, due to a linear algebra error.
While not definitive, this may be an indication of a poor mask.
Copy link
Member

Choose a reason for hiding this comment

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

Copy link
Member Author

Choose a reason for hiding this comment

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

We could do that, but why raise an error instead of issue a warning? The whole point of this is not to crash. People who care about CompCor can worry about why, but if they don't, the masks don't matter at all.

Please inspect the {} contours above to ensure that they are located
in the white matter/CSF.
</p>
'''.format(self._header, 'magenta' if self._header[0] == 'a' else 'blue')

with open(self._out_report, 'w') as fobj:
fobj.write(snippet)

class RobustTCompCor(nac.TCompCor):

class RobustACompCorInputSpec(RetryCompCorInputSpecMixin, nac.CompCorInputSpec):
pass


class RobustACompCorOutputSpec(reporting.ReportCapableOutputSpec, nac.CompCorOutputSpec):
pass


class RobustACompCor(RetryCompCorMixin, nac.ACompCor):
"""
Runs tCompCor several times if it suddenly fails with
https://github.com/poldracklab/fmriprep/issues/940
Runs aCompCor several times if it suddenly fails with
https://github.com/poldracklab/fmriprep/issues/776

Warns by default, rather than failing, on linear algebra errors.
https://github.com/poldracklab/fmriprep/issues/1433

"""
input_spec = RobustACompCorInputSpec
output_spec = RobustACompCorOutputSpec

def _run_interface(self, runtime):
failures = 0
while True:
try:
runtime = super(RobustTCompCor, self)._run_interface(runtime)
break
except LinAlgError:
failures += 1
if failures > 10:
raise
start = (failures - 1) * 10
sleep(randint(start + 4, start + 10))

return runtime
class RobustTCompCorInputSpec(RetryCompCorInputSpecMixin, nac.TCompCorInputSpec):
pass


class RobustTCompCorOutputSpec(reporting.ReportCapableOutputSpec, nac.TCompCorOutputSpec):
pass


class RobustTCompCor(RetryCompCorMixin, nac.TCompCor):
"""
Runs tCompCor several times if it suddenly fails with
https://github.com/poldracklab/fmriprep/issues/776

Warns by default, rather than failing, on linear algebra errors.
https://github.com/poldracklab/fmriprep/issues/1433

"""
input_spec = RobustTCompCorInputSpec
output_spec = RobustTCompCorOutputSpec
10 changes: 10 additions & 0 deletions fmriprep/viz/config.json
Original file line number Diff line number Diff line change
Expand Up @@ -98,6 +98,16 @@
"title": "ROIs in BOLD space",
"description": "Brain mask calculated on the BOLD signal (red contour), along with the masks used for a/tCompCor.<br />The aCompCor mask (magenta contour) is a conservative CSF and white-matter mask for extracting physiological and movement confounds. <br />The fCompCor mask (blue contour) contains the top 5% most variable voxels within a heavily-eroded brain-mask."
},
{
"name": "epi/warn_tcompcor",
"file_pattern": "func/.*_acompcor\\.",
"raw": true
},
{
"name": "epi/warn_tcompcor",
"file_pattern": "func/.*_tcompcor\\.",
"raw": true
},
{
"name": "epi_mean_t1_registration/flirt",
"file_pattern": "func/.*_flirtnobbr",
Expand Down
17 changes: 15 additions & 2 deletions fmriprep/workflows/bold/confounds.py
Original file line number Diff line number Diff line change
Expand Up @@ -180,12 +180,13 @@ def init_bold_confs_wf(mem_gb, metadata, name="bold_confs_wf"):
# a/t-CompCor
tcompcor = pe.Node(
TCompCor(components_file='tcompcor.tsv', header_prefix='t_comp_cor_', pre_filter='cosine',
save_pre_filter=True, percentile_threshold=.05),
save_pre_filter=True, percentile_threshold=.05, failure_mode='NaN',
generate_report=True),
name="tcompcor", mem_gb=mem_gb)

acompcor = pe.Node(
ACompCor(components_file='acompcor.tsv', header_prefix='a_comp_cor_', pre_filter='cosine',
save_pre_filter=True),
save_pre_filter=True, failure_mode='NaN', generate_report=True),
name="acompcor", mem_gb=mem_gb)

# Set TR if present
Expand Down Expand Up @@ -220,6 +221,16 @@ def init_bold_confs_wf(mem_gb, metadata, name="bold_confs_wf"):
name='ds_report_bold_rois', run_without_submitting=True,
mem_gb=DEFAULT_MEMORY_MIN_GB)

ds_report_warn_acompcor = pe.Node(
DerivativesDataSink(suffix='acompcor'),
name='ds_report_warn_acompcor', run_without_submitting=True,
mem_gb=DEFAULT_MEMORY_MIN_GB)

ds_report_warn_tcompcor = pe.Node(
DerivativesDataSink(suffix='tcompcor'),
name='ds_report_warn_tcompcor', run_without_submitting=True,
mem_gb=DEFAULT_MEMORY_MIN_GB)

def _pick_csf(files):
return files[0]

Expand Down Expand Up @@ -300,6 +311,8 @@ def _pick_wm(files):
(acc_msk, mrg_compcor, [('out', 'in2')]),
(mrg_compcor, rois_plot, [('out', 'in_rois')]),
(rois_plot, ds_report_bold_rois, [('out_report', 'in_file')]),
(acompcor, ds_report_warn_acompcor, [('out_report', 'in_file')]),
(tcompcor, ds_report_warn_tcompcor, [('out_report', 'in_file')]),
])

return workflow
Expand Down