[RTM] Add recon-all preprocessing #347
Conversation
fmriprep/workflows/anatomical.py
Outdated
| reconall = pe.Node( | ||
| freesurfer.ReconAll( | ||
| args='-parallel -openmp {:d}'.format(cpu_count())), | ||
| name='Reconstruction2') |
There was a problem hiding this comment.
Having an issue here:
traits.trait_errors.TraitError: The trait 'subjects_dir' of a ReconAllInputSpec instance is an existing
directory name, but the path '/scratch/workflow_enumerator/01/t1w_preprocessing/Reconstruction' does not exist.
But... it's still running Reconstruction2, and I'm getting outputs. Is this a nipype bug?
There was a problem hiding this comment.
Nope, actually it's a segfault.
There was a problem hiding this comment.
The erroring bit:
#-----------------------------------------
#@# Parcellation Stats 2 lh Thu Jan 26 17:41:39 UTC 2017
/scratch/workflow_enumerator/01/t1w_preprocessing/Reconstruction/recon_all/scripts
mris_anatomical_stats -th3 -mgz -cortex ../label/lh.cortex.label -f ../stats/lh.aparc.a2009s.stats -b -a ../label/lh.aparc.a2009s.annot -c ../label/aparc.annot.a2009s.ctab recon_all lh white
#-----------------------------------------
#@# Parcellation Stats 2 rh Thu Jan 26 17:41:39 UTC 2017
/scratch/workflow_enumerator/01/t1w_preprocessing/Reconstruction/recon_all/scripts
mris_anatomical_stats -th3 -mgz -cortex ../label/rh.cortex.label -f ../stats/rh.aparc.a2009s.stats -b -a ../label/rh.aparc.a2009s.annot -c ../label/aparc.annot.a2009s.ctab recon_all rh white
Waiting for PID 509 of (509 512) to complete...
/opt/freesurfer/bin/reconbatchjobs: line 81: 509 Segmentation fault (core dumped) exec $JOB >> $LOG 2>&1
Waiting for PID 512 of (509 512) to complete...
/opt/freesurfer/bin/reconbatchjobs: line 81: 512 Segmentation fault (core dumped) exec $JOB >> $LOG 2>&1
PIDs (509 512) completed and logs appended.
Linux 336cda4cfeae 4.4.0-59-generic #80-Ubuntu SMP Fri Jan 6 17:47:47 UTC 2017 x86_64 x86_64 x86_64 GNU/Linux
I've seen mris_anatomical_stats segfault issues before. I suspect it's because of the huge chunk of skull left by AFNI skull stripping (which is probably because of the low-res test data).
Short term: we'll kill the mris_anatomical_stats outputs, at least for testing. The outputs aren't used downstream in recon-all.
Medium term: I'll try to get a minimal reproduction (and better logs) for the FreeSurfer folks. Even if it fails, it shouldn't be segfaulting.
Somewhat relevant: There's no way to actually run these on CircleCI. I'll hook up the --no-freesurfer option for the tests, but we should figure out a way to test and validate FreeSurfer bits (at least) before releases.
fmriprep/workflows/anatomical.py
Outdated
| autorecon1 = pe.Node( | ||
| freesurfer.ReconAll( | ||
| directive='autorecon1', | ||
| args='-noskullstrip -parallel -openmp {:d}'.format(cpu_count())), |
There was a problem hiding this comment.
Please use the value from --nthreads
There was a problem hiding this comment.
Also set the number of expected threads for the scheduler (See https://github.com/poldracklab/fmriprep/blob/master/fmriprep/workflows/anatomical.py#L253)
effigies
force-pushed
the
recon_all
branch
6 times, most recently
from
9dd3826
to
a8cbcb7
Compare
|
Using https://openfmri.org/dataset/ds000005/ with Could someone have a look at this error? It produces the expected output file, but I can't figure out why it seems to think there should exist a Going to try running again in an absolutely clean scratch directory, but that's going to take a while, so if you see something right away, cool. |
|
Okay, so I'm having two major issues. The first is the above, where for whatever reason the reporting works fine, but throws an error. The second is I can't run Though that may be because ANTS skullstripping seems a bit over-zealous. |
effigies
force-pushed
the
recon_all
branch
6 times, most recently
from
0a920bc
to
40142ae
Compare
|
Report: sub-01.html BIDS directory after running: $ tree -L 4 /data/bids/ds000005_R2.0.1 | head -n 27
/data/bids/ds000005_R2.0.1
├── CHANGES
├── dataset_description.json
├── derivatives
│ └── freesurfer
│ └── sub-01
│ ├── label
│ ├── mri
│ ├── scripts
│ ├── stats
│ ├── surf
│ ├── tmp
│ ├── touch
│ └── trash
├── participants.tsv
├── README
├── sub-01
│ ├── anat
│ │ ├── sub-01_inplaneT2.nii.gz
│ │ └── sub-01_T1w.nii.gz
│ └── func
│ ├── sub-01_task-mixedgamblestask_run-01_bold.nii.gz
│ ├── sub-01_task-mixedgamblestask_run-01_events.tsv
│ ├── sub-01_task-mixedgamblestask_run-02_bold.nii.gz
│ ├── sub-01_task-mixedgamblestask_run-02_events.tsv
│ ├── sub-01_task-mixedgamblestask_run-03_bold.nii.gz
│ └── sub-01_task-mixedgamblestask_run-03_events.tsv
... |
|
At this point, we have reconstruction, saving derivatives, and reporting. At present, I'm accepting whatever external skullstripping is used, ANTs or AFNI. If we want to only swap out the FreeSurfer skullstripping if ANTs is specified, we could do that. The following features are coded but not tested:
Talking with @chrisfilo, perhaps the T2 stuff should be stripped out for now, and another PR can pull them in when we can get some proper testing. The same could go for multiple T1s or The whole thing depends on nipy/nipype#1790, which serves as a decent test of that PR. |
|
Yes - lets keep T2s for another PR. We should get some test data soon.
…On Tue, Feb 7, 2017 at 2:22 PM, Chris Markiewicz ***@***.***> wrote:
At this point, we have reconstruction, saving derivatives, and reporting.
At present, I'm accepting whatever external skullstripping is used, ANTs
or AFNI. If we want to only swap out the FreeSurfer skullstripping if ANTs
is specified, we could do that.
The following features are coded but not tested:
- Multiple T1w images (haven't yet dug up a dataset to test this on)
- T2 images (no <1.2mm T2s to test on)
- Also, at present accepting inplaneT2 (because those were
available in the above dataset), but I think we should actually be looking
for images labeled T2w?
- -hires (Waiting in a Sherlock queue to run on BEAST data)
Talking with @chrisfilo <https://github.com/chrisfilo>, perhaps the T2
stuff should be stripped out for now, and another PR can pull them in when
we can get some proper testing. The same could go for multiple T1s or
-hires, if we just want to get basic support in as quickly as possible.
The whole thing depends on nipy/nipype#1790
<nipy/nipype#1790>, which serves as a decent test
of that PR.
—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
<#347 (comment)>,
or mute the thread
<https://github.com/notifications/unsubscribe-auth/AAOkp2Aoa5233OGZdCoyZcfXB1yptqxhks5raO67gaJpZM4LuLEI>
.
|
|
Are you planning to map the EPI timeseries to fsaverage in a separate PR? |
|
Yes, since it's a separate process to |
effigies
force-pushed
the
recon_all
branch
3 times, most recently
from
7c4179f
to
eb37f93
Compare
| def detect_inputs(t1w_list, default_flags=''): | ||
| from nipype.utils.filemanip import filename_to_list | ||
| import nibabel as nib | ||
| t1w_list = filename_to_list(t1w_list) |
There was a problem hiding this comment.
If t1w is undefined, this function will have issues. However since it's a mandatory input to the workflow, it also should cause general failure, so I haven't specifically handled that case.
fmriprep/interfaces/bids.py
Outdated
| @@ -23,6 +23,11 @@ | |||
| from fmriprep.utils.misc import collect_bids_data, make_folder | |||
|
|
|||
| LOGGER = logging.getLogger('interface') | |||
| BIDS_NAME = re.compile( | |||
| '^(.*\/)?(?P<subject_id>sub-[a-zA-Z0-9]+)(_(?P<ses_id>ses-[a-zA-Z0-9]+))?' | |||
There was a problem hiding this comment.
Note that this changes from only matching sub-... at the start of string to also matching after a /.
| autorecon1_flags.append('-hires') | ||
| reconall_flags.append('-hires') | ||
| return (t1w_outs, ' '.join(autorecon1_flags), | ||
| ' '.join(reconall_flags)) |
There was a problem hiding this comment.
This approach doesn't permit type checking the flags that are passed. -noskullstrip and -hires aren't ReconAll options, but this will eventually also produce -T2 and -T2pial options, which are supported directly by nipype.
Not sure if this will bother anybody.
There was a problem hiding this comment.
Doesn't -hires need to be entered directly before a corresponding -i with a T1 volume for it to work?
There was a problem hiding this comment.
Not that I can tell from reading the Csh source, it looks like it's set globally. Do you have a source suggesting this?
There was a problem hiding this comment.
Though, it does look like it should come with an experts file, which we should probably produce. I guess just go with their mris_inflate -n 15?
There was a problem hiding this comment.
Did you check in the logs that the right resolution is used internally?
BTW it actually runs without an experts file (at least in my experience).
There was a problem hiding this comment.
Yeah, it does run without the experts file.
Not sure what to look for in the log, but (from a currently running test):
>>> nib.load('/local-scratch/cjmarkie/12650555/workflow_enumerator/sub-387/t1w_preprocessing/Reconstruction/sub-387/mri/T1.mgz').header.get_zooms()
(0.89840001, 0.89840001, 0.89840001)There was a problem hiding this comment.
Pinged the FreeSurfer list about -hires best practices. My guess is that 15 max iterations should work for anything -- if something works well with a max of 10, it should stop by 10, but we should give smaller voxels a chance. 760d9de has the update I'd make if we want to go with that strategy. That seems like it's probably a good idea, since we'll want to support people with smaller voxels than we currently have access to.
| match = BIDS_NAME.search(in_file) | ||
| params = match.groupdict() if match is not None else {} | ||
| return tuple(map(params.get, ['subject_id', 'ses_id', 'task_id', | ||
| 'acq_id', 'rec_id', 'run_id'])) |
There was a problem hiding this comment.
Would something like this be more generally useful? I could put it in a better home if we might want to reuse.
There was a problem hiding this comment.
This could be useful for pybids. I've opened an issue bids-standard/pybids#31
fmriprep/workflows/anatomical.py
Outdated
| autorecon1.interface._can_resume = False | ||
| autorecon1.interface.num_threads = nthreads | ||
| if 'FREESURFER_SUBJECTS' in os.environ: | ||
| autorecon1.inputs.subjects_dir = os.getenv('FREESURFER_SUBJECTS') |
There was a problem hiding this comment.
This doesn't get used in my docker runs; do we even care about supporting existing FREESURFER_SUBJECTS directories, or is the derivatives folder the correct BIDS-y way to do this? A command-line option is always a possibility, as well.
fmriprep/workflows/anatomical.py
Outdated
| mask = nib.load(skullstripped) | ||
| bmask = new_img_like(mask, mask.get_data() > 0) | ||
| resampled_mask = resample_to_img(bmask, img, 'nearest') | ||
| masked_image = new_img_like(img, img.get_data() * resampled_mask.get_data()) |
There was a problem hiding this comment.
Mentioned this to Chris, but for anybody else who's reviewing, just note that the choice of binarizing before resampling will be slightly more prone to declaring voxels "not brain" than after. This is because, at the boundary, an interpolated value will always be > 0, while a nearest neighbor in a resampled binary map may be 0.
Just in case anybody wants to object.
fmriprep/workflows/anatomical.py
Outdated
| masked_image = new_img_like(img, img.get_data() * resampled_mask.get_data()) | ||
| masked_image.to_filename(bm_auto) | ||
|
|
||
| copyfile(bm_auto, bm, copy=True, use_hardlink=True) |
There was a problem hiding this comment.
Hardlinks should be fine, right? Seems a waste to copy.
There was a problem hiding this comment.
Brain masks are so small and hard links are not supported on all systems. Will this throw an exception if hard links are not supported?
There was a problem hiding this comment.
Yeah, the fallback behavior would be to copy.
fmriprep/workflows/anatomical.py
Outdated
| input_names=['dirname', 'basename'], | ||
| output_names=['path']), | ||
| name='SubjectDir', | ||
| run_without_submitting=True) |
There was a problem hiding this comment.
Is there a better way to get the full subject directory? It's 11 lines for os.path.join.
fmriprep/workflows/anatomical.py
Outdated
| base_directory=settings['output_dir'], | ||
| container='derivatives'), | ||
| name='FSSubjOut', | ||
| run_without_submitting=True) |
There was a problem hiding this comment.
Used DataSink instead of DerivativesDataSink because the latter doesn't like directories.
Disable version check in info.py for now
|
Derivatives properly copying. Rebased against Any more comments? |
fmriprep/workflows/anatomical.py
Outdated
| run_without_submitting=True) | ||
| # Copy fsaverage subject into derivatives/freesurfer | ||
| setattr(subject_out.inputs, 'freesurfer.@fsaverage', | ||
| '/opt/freesurfer/subjects/fsaverage') |
There was a problem hiding this comment.
Working on bbregister/mri_vol2surf, I'm realizing we don't want this inside t1w_preprocessing -- or really anything below workflow_enumerator -- as then it would be re-run for every subject.
Does it make sense to make something like:
class BIDSFreeSurferDir(BaseInterface):
inputs = (freesurfer_home, derivatives_dir)
outputs = subjects_dir
def _run_interface(self, runtime):
# Create subjects directory
# Copy fsaverageAnd then pass subjects_dir through to workflows that need access to it? We would then be recon-all'ing directly in derivatives/freesurfer/sub-{subj_id}/, which seems okay to me.
WDYT, @chrisfilo? This can wait, but I think we'll have to deal with it in the next PR (where bbregister will need to know about the common subjects directory to do the registration) if we don't, now.
effigies
force-pushed
the
recon_all
branch
3 times, most recently
from
da03c28
to
127651b
Compare
|
Okay, so I've tested the following cases on Sherlock:
Pushed a commit to address this last case and not overwrite brainmasks. I believe we agreed that existing recon directories should be treated as "correct", and if someone wants to start over, they can delete the recon directory. |
|
+1
…On Thu, Feb 16, 2017 at 1:04 PM, Chris Markiewicz ***@***.***> wrote:
Okay, so I've tested the following cases on Sherlock:
1. New derivatives directory - copies fsaverage, creates new subject
2. Existing fsaverage, new subject - leaves fsaverage untouched,
creates new subject
3. Existing fsaverage and subject - leaves fsaverage untouched, runs
ReconAll without error (also not doing anything)
4. Existing fsaverage and subject, deleted nipype scratch directory -
leaves fsaverage untouched, runs ReconAll without error (overwriting
brainmasks)
Pushed a commit to address this last case and not overwrite brainmasks. I
believe we agreed that existing recon directories should be treated as
"correct", and if someone wants to start over, they can delete the recon
directory.
—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
<#347 (comment)>,
or mute the thread
<https://github.com/notifications/unsubscribe-auth/AAOkp1JojssU8zagT2e3RPzCB2YRPpyrks5rdLn0gaJpZM4LuLEI>
.
|
|
Alright. No errors, and existing brainmasks were left alone in case 4. This is ready for final review and merge. |
|
Excellent work! |
The goal of this PR is to get recon-all capabilities into FMRIPREP. Using its outputs in other steps in the pipeline will be addressed in other PRs.
Goals from #101:
-hiresfor T1 voxel size < 1mmmri_converting ANTS output tomri/brainmask(.auto).mgzsqrt(T1*T2)bias field correction#101 issues for a separate PR (for easy reference later):
bbregisterto register functionals to T1