ENH: Temporary integration of FSL eddy #144
Conversation
Dockerfile
Outdated
| ENV FSLDIR="/opt/fsl-6.0.1" \ | ||
| PATH="/opt/fsl-6.0.1/bin:$PATH" \ | ||
| FSLOUTPUTTYPE="NIFTI_GZ" | ||
|
|
||
| RUN apt-get update -qq && \ | ||
| apt-get install -y -q --no-install-recommends \ | ||
| bc \ | ||
| dc \ | ||
| file \ | ||
| libfontconfig1 \ | ||
| libfreetype6 \ | ||
| libgl1-mesa-dev \ | ||
| libglu1-mesa-dev \ | ||
| libgomp1 \ | ||
| libice6 \ | ||
| libxcursor1 \ | ||
| libxft2 \ | ||
| libxinerama1 \ | ||
| libxrandr2 \ | ||
| libxrender1 \ | ||
| libxt6 \ | ||
| python \ | ||
| wget && \ | ||
| apt-get clean && \ | ||
| rm -rf /var/lib/apt/lists/* && \ | ||
| echo "Downloading FSL ..." && \ | ||
| wget -q http://fsl.fmrib.ox.ac.uk/fsldownloads/fslinstaller.py && \ | ||
| chmod 775 fslinstaller.py && \ | ||
| /usr/bin/python fslinstaller.py -d /opt/fsl-6.0.1 -V 6.0.1 -q | ||
|
|
There was a problem hiding this comment.
Is it necessary to upgrade to FSL 6? Unless the interface has changed dramatically, I would leave that for a different PR
There was a problem hiding this comment.
I put this in because the current docker container didn't have eddy_openmp which the node is looking for (in terms of that specific command and also parameters). I thought just installing FSL 6 is the best way to go for this, I'm open to changing it to something else that works though.
There was a problem hiding this comment.
The Nipype interface for EddyCorrect (the older command for eddy) might work here instead. eddy_openmp introduced better susceptibility distortion correction but that wouldn't be necessary if we use sdcflows for that.
There was a problem hiding this comment.
Thanks, in that case I can just use that. EddyCorrect doesn't need the index file and acqp inputs, so would it be best to scrap those nodes for now and just use EddyCorrect?
There was a problem hiding this comment.
One thing I'm unsure of is whether the old eddy_correct can handle multi-shell or high b-value data. But our test case (the travelling human phantoms) is single-shell.
There was a problem hiding this comment.
I edited it to just use eddy_correct for now, but the other functions needed for eddy_openmp are still there in case we want to upgrade. Let me know what you guys think.
dmriprep/workflows/dwi/util.py
Outdated
| if metadata.get("PhaseEncodingDirection"): | ||
| pe_dir = metadata.get("PhaseEncodingDirection") | ||
| else: | ||
| pe_dir = metadata.get("PhaseEncodingAxis") |
There was a problem hiding this comment.
Is this valid?
(still think this should be passed as inputs)
There was a problem hiding this comment.
This was an old backup check that I don't believe the PhaseEncodingAxis part is accurate now that you mention it. I do still find getting the PhaseEncodingDirection from each subject's metadata useful if its not manually passed as an input. So is it alright to change to use inputs, otherwise revert to PhaseEncodingDirection?
dmriprep/workflows/dwi/util.py
Outdated
| else: | ||
| total_readout = total_readout_time | ||
|
|
||
| # Construct the acqp file lines |
There was a problem hiding this comment.
lines? AFAICT, this will generate just one line right?
There was a problem hiding this comment.
Right, double checked and it generates one line.
There was a problem hiding this comment.
sdcflows also has a function for getting/calculating this: https://github.com/nipreps/sdcflows/blob/66b3c3009b9d951b19f36902c51e379d36e292f2/sdcflows/utils/epimanip.py#L14
There was a problem hiding this comment.
Thanks, added sdcflows' function in e3d347c
|
Hello @slimnsour, Thank you for updating! Cheers! There are no style issues detected in this Pull Request. 🍻 To test for issues locally, Comment last updated at 2021-02-06 21:20:13 UTC |
dmriprep/workflows/dwi/util.py
Outdated
| name="gen_idx", | ||
| ) | ||
|
|
||
| eddy_quad = pe.Node(EddyQuad(verbose=True), name="eddy_quad") |
There was a problem hiding this comment.
We won't need to include eddy_quad. That's used for creating QC pages but this pipeline will have its own version.
There was a problem hiding this comment.
Thanks for the catch, fixed in 6605fe1
|
Okay, I've sent some more concrete suggestions to slimnsour#1 - please review and merge if you agree. EDIT |
|
BTW - we will need to also tinker with the parameters so that we can run a fast (and inaccurate) eddy process on Circle (i.e., with the |
|
Sorry for the delay on this. I figured out that the issue with the tests was that there were outdated python2 print statements in code used by And as for the |
|
I believe the current build fails are related to current issues with svgutils' new version (see nipreps/niworkflows#595) messing with niworkflows, it has the same |
|
Yup, can you check whether the version of svgutils being installed on the docker build step of CircleCI has changed here w.r.t. the latest successful master? |
|
From sshing into the circleci build I can confirm that it has the problematic svgutils, 0.3.2. Also, going back to 0.3.1 no longer causes the int64 issues. Would the best move here be blocking 0.3.2 from being installed until a fix is put in? |
|
If you want to be a good citizen, you probably want to report this in their repo, pin svgutils<0.3.2, and then undo the pinning if it's not a bug and rather a misuse no longer allowed after 0.3.1. |
|
Thanks for mentioning the issue. It seems it was fixed in btel/svg_utils#57 (comment) (which I have also tested and can confirm it works) but a new release isn't out yet. Either way I'll exclude the problematic release. |
|
Thanks @slimnsour! I just checked the QC reports. ds000206 looks good but the mask in ds001771 looks wonky. Should we try rerunning that workflow? |
|
Good call, not sure what happened there. I'll try rerunning. Edit: I just checked my own test output and it seemed fine too. |
|
Rerun finished and mask looks good! Let me know if there's anything else to do, otherwise I think we can merge this. |
dmriprep/workflows/dwi/eddy.py
Outdated
| included in FSL {EddyCorrect().version} [@eddy]. | ||
| """ | ||
|
|
||
| eddy_correct = pe.Node(EddyCorrect(), name="eddy_correct",) |
There was a problem hiding this comment.
I've been checking and EddyCorrect only returns the data after correction, which is a big issue. We will need the head-motion parameters, at the very least, for further processing.
Sorry to push back on this. I went to our earlier conversation to check why we pivoted that way. We (or I) thought that FSL 5 does not have eddy_openmp - so fell back to the deprecated version.
However:
/opt/fsl-5.0/bin/eddy_openmp
***************************************************
The following COMPULSORY options have not been set:
--imain File containing all the images to estimate distortions for
--mask Mask to indicate brain
--index File containing indices for all volumes in --imain into --acqp and --topup
--acqp File containing acquisition parameters
--bvecs File containing the b-vectors for all volumes in --imain
--bvals File containing the b-values for all volumes in --imain
--out Basename for output
***************************************************
Part of FSL (ID: 5.0.11)
At least 5.0.11 has eddy_openmp.
There was a problem hiding this comment.
My mistake. Upon opening the code for eddy_correct, it looks like it was just a wrapper for flirt. Thankfully, @slimnsour has kept the functions for writing the index and acquisition parameters file. v5.0.11 introduced the ability to do slice-to-volume motion correction. This could be useful if the dwi json includes slice timing values.
There was a problem hiding this comment.
Sorry for the delay on this. I figured out that the issue with the tests was that there were outdated python2 print statements in code used by
eddy_correct. I created an updated file that letseddy_correctrun correctly and it is used to overwrite the base FSL one through our.dockerfolder.And as for the
--sloppyimplementation, I'm not sure what else can be done with the limitededdy_correctparameters to make it faster. I'd be happy to implement any solutions if you guys have any (maybe a parameter that skipseddy_correctaltogether?)
If migrating to eddy_openmp, one suggestion could be decreasing the number of iterations with niter.
There was a problem hiding this comment.
Thanks for the feedback guys, attempting to implement this now with 60efdaf.
It seems FSL 5.0.11 is not available using apt-get, so I am putting in an installation step in the Dockerfile to get eddy_openmp. And I set the --sloppy tag to make niter only 1 (which may be overkill but should be much faster).
There was a problem hiding this comment.
An alternative to installing FSL 5.0.11 could be to just include the eddy_openmp script in the .docker folder (similar to how topup was added). But I'm not sure if there are added dependencies that are also required.
There was a problem hiding this comment.
(or does eddy_openmp use it internally?)
There was a problem hiding this comment.
Yes, the masking workflow still uses BET. I don't think eddy_openmp uses BET internally. It would rely on the mask that we pass it.
There was a problem hiding this comment.
Good suggestion, I've tested just using the one script and it seems to be working just fine so it's pushed in f463e66. I was worried that eddy_openmp would have dependencies that we couldn't individually account for so I initially just installed the entire FSL 5.0.11.
dmriprep/workflows/dwi/eddy.py
Outdated
| included in FSL {Eddy().version} [@eddy]. | ||
| """ | ||
| eddy = pe.Node( | ||
| Eddy(repol=True, cnr_maps=True, residuals=True, method="jac"), |
There was a problem hiding this comment.
| Eddy(repol=True, cnr_maps=True, residuals=True, method="jac"), | |
| Eddy(), |
Not necessary to output cnr maps or residual image. method is also jacobians by default. Not sure how we want to handle outlier slices/volumes. This would decide whether we should keep repol=True or not. It could definitely be False in sloppy mode.
|
It seems the ds000206 dataset doesn't have the necessary info for running Do you guys have any suggestions on how we can go about this? Perhaps a parameter that allows for a default value to be used if no trt is calculated would work here. For reference, the |
|
I wonder we could drop a dummy value and raise a warning in these cases. Reading the "Does it matter what I put into my --acqp file?" section here, it seems like it might work in many cases. Could you give it a try with this dataset? Maybe we can look at the outputs to see whether it does something crazy? |
|
Thanks for the suggestion @arokem, I think that's a fair addition. I have implemented this now with a backup value of |
|
If the parameter has no influence over the outcome, then I think dropping a dummy value is fine. However, if it does, then I would fix our test dataset, write the dummy value there, and let dMRIPrep break when people do not specify their dataset in full. That will force them to write the dummy value themselves, which hopefully will bring them to consider two things:
Let's take advantage of BIDS here! |
| @@ -298,7 +298,7 @@ jobs: | |||
| - /tmp/ds000206/work | |||
| - run: | |||
| name: Run full diffusion workflow on ds000206 | |||
| no_output_timeout: 2h | |||
| no_output_timeout: 4h | |||
There was a problem hiding this comment.
Is there a way of running this in some sort of "sloppy/fast" mode? 4h is a lot!
There was a problem hiding this comment.
4h is a lot for future testing, and this was already cut down by decreasing the --niter tag from 5 to 1 in sloppy mode as per @josephmje's suggestion. I will try to decrease other parameters like --repol (and perhaps --mbs_niter?).
There was a problem hiding this comment.
I guess :D - I have no idea of how to control eddy.
dmriprep/workflows/dwi/base.py
Outdated
| ) | ||
|
|
||
| eddy_report = pe.Node( | ||
| SimpleBeforeAfter(before_label="Distorted", after_label="Eddy Corrected",), |
There was a problem hiding this comment.
I wonder whether there is a better reportlet. Although you can see that the average image seems to have a better SNR after eddy, do you think that some sort of standard deviation map across orientations or similar will give us a better picture of what this step did (and whether it did it well)?
WDTY @arokem ?
There was a problem hiding this comment.
The developers of eddy seem to like showing an animation of flipping through all the different directions rather quickly.
There was a problem hiding this comment.
Another option is to show before and after mean maps across all directions (is that what this does?).
Unfortunately, I think that an std map would capture both artifactual and real variance (e.g., due to diffusion anisotropy).
There was a problem hiding this comment.
So, now that you made me check what this exactly does I can tell that it shows the b=0 reference we generate (averaging all b=0 in the dataset after regressing out signal drift) and the first volume of the eddy-corrected output (which happens to be a b=0, as it is the most common case). So at the very least, we want to plot two b=0 references generated the same way.
I'm not sure that we are going to see much averaging all directions... however, I believe the std map would be still useful - since you have the reference std map (that would be the "Before" on the flickering mosaic) and then you plot the std map after correction (the "After"), adjusted to have consistent intensity ranges, we should see how some of the variability (that from the artifact) goes away. Because of the structure of the artifact, I would expect better SNR (less noise in the "after" plot) at the edges of the WM, esp. around the ventricles. WDYT?
There was a problem hiding this comment.
For the averaging, I imagine the "pre" image would have a lot more blur than the "post"?
But for the std image: I think that you have me convinced.
There was a problem hiding this comment.
Sounds like a good idea! I'd just need some direction so I can implement it - where is the std map input coming from in the workflow?
There was a problem hiding this comment.
This paper expands on what FSL currently outputs for QCing eddy: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6264528/pdf/main.pdf
There was a problem hiding this comment.
Just checking up on this @oesteban, you mentioned in the meeting that you would be able to throw something together for the reports you guys were interested in?
There was a problem hiding this comment.
Yup, sorry I've been booked up. Let's merge this and address the reportlet on a different PR.
|
Please rebase on upstream/master when #151 is merged :) |
dmriprep/config/reports-spec.yml
Outdated
| @@ -91,6 +91,13 @@ sections: | |||
| the final images are resampled using cubic B-Spline interpolation. | |||
| static: false | |||
| subtitle: Alignment of diffusion and anatomical MRI data (surface driven) | |||
| - bids: {datatype: figures, desc: eddy, suffix: dwi} | |||
| caption: Head-motion and eddy current corrected diffusion data | |||
| description: Geometrical distortions derived from the so-called Eddy-currents, | |||
There was a problem hiding this comment.
| description: Geometrical distortions derived from the so-called Eddy-currents, | |
| description: Eddy current-induced geometric distortions |
There was a problem hiding this comment.
@slimnsour would you accept this suggestion?
There was a problem hiding this comment.
Yup, fixed in caba8c3
dmriprep/config/reports-spec.yml
Outdated
| @@ -91,6 +91,13 @@ sections: | |||
| the final images are resampled using cubic B-Spline interpolation. | |||
| static: false | |||
| subtitle: Alignment of diffusion and anatomical MRI data (surface driven) | |||
| - bids: {datatype: figures, desc: eddy, suffix: dwi} | |||
There was a problem hiding this comment.
This section should be moved before the alignment with the anatomical data.
There was a problem hiding this comment.
Thanks, fixed in caba8c3
Just chiming in to confirm that eddy does not use total readout time if we don't provide the fieldmap or topup arguments. I've tested this out with a few different dummy values in the past and would get the same result. We may face an issue when attempting to apply a true fieldmap to the dwi data and total readout time can't be calculated. But this can just be resolved by sending a warning/erroring out. |
dmriprep/config/reports-spec.yml
Outdated
| @@ -91,6 +91,13 @@ sections: | |||
| the final images are resampled using cubic B-Spline interpolation. | |||
| static: false | |||
| subtitle: Alignment of diffusion and anatomical MRI data (surface driven) | |||
| - bids: {datatype: figures, desc: eddy, suffix: dwi} | |||
Co-authored-by: Oscar Esteban <code@oscaresteban.es>
Co-authored-by: Oscar Esteban <code@oscaresteban.es>
Co-authored-by: Michael Joseph <josephmje.22@gmail.com>
Porting over from our old dmriprep project. I also had to install a full version of FSL in the Dockerfile in order for the pipeline to have access to eddy.