FSL TOPUP changes behavior at nipype 1.0.3

[billbrod]

Summary

The way the interface to FSL TOPUP works changes at nipype 1.0.3. Previously, the defualt behavior was not to pass a flag for the fwhm argument, which meant that TOPUP used its default argument (8, 6, 4, 3, 3, 2, 1, 0, 0) for fwhm. Note that is not what the FSL documentation states is the default behavior (the documentation says the default is to use a value of 8), but checking the log file produced by FSL shows that this is what happens. (See this email to the FSL listserv for some more details).

After version 1.0.3, nipype explicitly passes --fwhm 8 as its default behavior, conforming with what FSL says in the docs, but not conforming with TOPUP's actual default behavior. This can be changed to a single number, but cannot be changed to a list (related to issue #908) and so the previous behavior (passing --fwhm 8, 6, 4, 3, 3, 2, 1, 0, 0) is impossible to achieve.

I don't understand why this is the default FSL behavior, but from a little bit of playing around, it does appear to result in less distorted images (compared to --fwhm 0) which are unblurred (unlike --fwhm 8).

We noticed this by comparing the outputs of nipype 0.13 and 1.1.0. I'm saying the change happens at 1.0.3 based on the email to the FSL listserv linked above. This is all using python 3.

Actual (nipype >= 1.0.3) behavior

command.txt produced by nipype

topup --config=b02b0.cnf --datain=/Volumes/server/Projects/BAIR/Analyses/visual/sub-som682/ses-nyu3t01/workflow/unwarp_distort/distortion_merged_encfile.txt --fwhm=8.000000 --imain=/Volumes/server/Projects/BAIR/Analyses/visual/sub-som682/ses-nyu3t01/workflow/merge_distort/distortion_merged.nii.gz --miter=5 --out=distortion_merged_base --iout=distortion_merged_corrected.nii.gz --fout=distortion_merged_field.nii.gz --jacout=jac --logout=distortion_merged_topup.log --rbmout=xfm --dfout=warpfield --miter=1 --splineorder=3 --subsamp=1 --warpres=10.000000

log produced by FSL

# 	name of 4D file with images
--imain=/Volumes/server/Projects/BAIR/Analyses/visual/sub-som682/ses-nyu3t01/workflow/merge_distort/distortion_merged.nii.gz
# name of text file with PE directions/times
--datain=/Volumes/server/Projects/BAIR/Analyses/visual/sub-som682/ses-nyu3t01/workflow/unwarp_distort/distortion_merged_encfile.txt
# 	base-name of output files (spline coefficients (Hz) and movement parameters)
--out=distortion_merged_base
# 	name of image file with field (Hz)
--fout=distortion_merged_field.nii.gz
# 	name of 4D image file with unwarped images
--iout=distortion_merged_corrected.nii.gz
# (approximate) resolution (in mm) of warp basis for the different sub-sampling levels, default 10
--warpres=10
# sub-sampling scheme, default 1
--subsamp=1
# 	FWHM (in mm) of gaussian smoothing kernel, default 8
--fwhm=8
# 	Max # of non-linear iterations, default 5
--miter=1
# Weight of regularisation, default depending on --ssqlambda and --regmod switches. See user documetation.
--lambda=0.005,0.001,0.0001,1.5e-05,5e-06,5e-07,5e-08,5e-10,1e-11
# If set (=1), lambda is weighted by current ssq, default 1
--ssqlambda=1
# Estimate movements if set, default 1 (true)
--estmov=1,1,1,1,1,0,0,0,0
# Minimisation method 0=Levenberg-Marquardt, 1=Scaled Conjugate Gradient, default 0 (LM)
--minmet=0,0,0,0,0,1,1,1,1
# Model for regularisation of warp-field [membrane_energy bending_energy], default bending_energy
--regmod=bending_energy
# Order of spline, 2->Qadratic spline, 3->Cubic spline. Default=3
--splineorder=3
# Precision for representing Hessian, double or float. Default double
--numprec=double
# Image interpolation model, linear or spline. Default spline
--interp=spline
# 	If set (=1), the images are individually scaled to a common mean, default 0 (false)
--scale=1
# 	If set (=1), the calculations are done in a different grid, default 1 (true)
--regrid=1

Expected behavior

We took the command from above and just removed the --fwhm 8 flag:

topup --config=b02b0.cnf --datain=/Volumes/server/Projects/BAIR/Analyses/visual/sub-som682/ses-nyu3t01/workflow/unwarp_distort/distortion_merged_encfile.txt --imain=/Volumes/server/Projects/BAIR/Analyses/visual/sub-som682/ses-nyu3t01/workflow/merge_distort/distortion_merged.nii.gz --miter=5 --out=distortion_merged_base --iout=distortion_merged_corrected.nii.gz --fout=distortion_merged_field.nii.gz --jacout=jac --logout=distortion_merged_topup.log --rbmout=xfm --dfout=warpfield --miter=1 --splineorder=3 --subsamp=1 --warpres=10.000000

log file produced by FSL

# 	name of 4D file with images
--imain=/Volumes/server/Projects/BAIR/Analyses/visual/sub-som682/ses-nyu3t01/workflow/merge_distort/distortion_merged.nii.gz
# name of text file with PE directions/times
--datain=/Volumes/server/Projects/BAIR/Analyses/visual/sub-som682/ses-nyu3t01/workflow/unwarp_distort/distortion_merged_encfile.txt
# 	base-name of output files (spline coefficients (Hz) and movement parameters)
--out=distortion_merged_base
# 	name of image file with field (Hz)
--fout=distortion_merged_field.nii.gz
# 	name of 4D image file with unwarped images
--iout=distortion_merged_corrected.nii.gz
# (approximate) resolution (in mm) of warp basis for the different sub-sampling levels, default 10
--warpres=10
# sub-sampling scheme, default 1
--subsamp=1
# 	FWHM (in mm) of gaussian smoothing kernel, default 8
--fwhm=8,6,4,3,3,2,1,0,0
# 	Max # of non-linear iterations, default 5
--miter=1
# Weight of regularisation, default depending on --ssqlambda and --regmod switches. See user documetation.
--lambda=0.005,0.001,0.0001,1.5e-05,5e-06,5e-07,5e-08,5e-10,1e-11
# If set (=1), lambda is weighted by current ssq, default 1
--ssqlambda=1
# Estimate movements if set, default 1 (true)
--estmov=1,1,1,1,1,0,0,0,0
# Minimisation method 0=Levenberg-Marquardt, 1=Scaled Conjugate Gradient, default 0 (LM)
--minmet=0,0,0,0,0,1,1,1,1
# Model for regularisation of warp-field [membrane_energy bending_energy], default bending_energy
--regmod=bending_energy
# Order of spline, 2->Qadratic spline, 3->Cubic spline. Default=3
--splineorder=3
# Precision for representing Hessian, double or float. Default double
--numprec=double
# Image interpolation model, linear or spline. Default spline
--interp=spline
# 	If set (=1), the images are individually scaled to a common mean, default 0 (false)
--scale=1
# 	If set (=1), the calculations are done in a different grid, default 1 (true)
--regrid=1

[djarecka]

@billbrod - thank you for reporting, the changes has been also reported in the comments of the PR that changed default values of topup.
I opened a new PR to fix it. You can check it out and let me know if you have suggestions.

[billbrod]

Thanks for your quick response and sorry I didn't see the existing PR. I'm not sure the best practice is in this situation -- should I close this issue and take any discussion over to the PR or continue the discussion here?

Your PR looks like it fixes the problem I mentioned, but I think we would also want to allow multiple values for the different TOPUP arguments (that is, implement the feature described in #908). I find it very confusing that FSL's documentation of its default behavior does not match its actual default behavior. I'm not sure whether it's better to let FSL handle the defaults or make nipype's default behavior match up with FSL's (the second has the advantage of making things more explicit and resilient to any changes in FSL's default behavior), but I think the documentation should be updated to accurately reflect what is happening. So, for instance, the description of the fwhm parameter shouldn't be "FWHM (in mm) of gaussian smoothing kernel (default in FSL 8)", but something along the lines of "FWHM (in mm) of gaussian smoothing kernel (default in FSL 8, 6, 4, 3, 3, 2, 1, 0, 0)". What do you think?

[djarecka]

@billbrod - thanks for additional comments. For now, I simply concentrated on returning to the previous behavior, it should be part of the new release on Mon. It wasn't sure how to provide a meaningful information, for example FWHM is an optional argument, so by default nothing is set, etc. Just decided that it's better to motivate people to go to FSL directly than provide misleading info. But if you have suggestions, please feel free to submit a PR.

Regarding #908, I don't think anyone is working on this right now (@oesteban, please correct me if I'm wrong). Would you like to give it a try?

[oesteban]

Yes, I think I never pushed any changes to #908.

[mgxd]

closed in #2637

[billbrod]

@djarecka Alright thanks, glad to see this change implemented! I'll look into #908 and see if I can contribute to it.

[djarecka]

sure, it will be in today release. And please feel free to ask any question if you want to try #908.