nipreps / fmriprep

fMRIPrep is a robust and easy-to-use pipeline for preprocessing of diverse fMRI data. The transparent workflow dispenses of manual intervention, thereby ensuring the reproducibility of the results.
https://fmriprep.org
Apache License 2.0
634 stars 293 forks source link

TOPUP issues in functionals, even though the topup field looks reasonable #2835

Closed cmpetty closed 11 months ago

cmpetty commented 2 years ago

We have some pretty severely distorted images because of a protocol mistake, but we're trying to make the best of them.

PEPolar images below: Screen Shot 2022-08-15 at 6 11 25 PM

The actual SDC correction looks reasonable on the pepolar images:

Screen Shot 2022-08-18 at 9 17 47 AM

However, when the field is applied to the functionals things go wrong:

Screen Shot 2022-08-18 at 9 18 26 AM

I've verified the PE directions and the readouts are correct. I've limited the number of topup volumes to use, just as a test, which didn't matter. I've tracked that the images look normal up until the "unwarp_wf" stage.

I'm using v21.0.2 via singularity containers with the command below: ${PACKDIR}/singularity/bin/singularity run --cleanenv -B $PACKDIR:/mnt/packages \ ${PACKDIR}/singularity/images/fmriprep-v21.0.2.sif \ /mnt/munin/Neacsiu/${EXPERIMENT}/Data/BIDS/TMSCR \ /mnt/munin/Neacsiu/${EXPERIMENT}/Data/BIDS/fmriprep \ participant --participant-label ${SUBJ} \ --nthreads 12 --omp-nthreads 12 --skip_bids_validation \ --use-aroma --aroma-melodic-dimensionality 30 \ --fs-no-reconall \ --output-spaces MNI152NLin2009cAsym:res-02 \ --work-dir /mnt/munin/Neacsiu/${EXPERIMENT}/Analysis/work \ --topup-max-vols 1 \

Outside of fmriprep i have run topup myself, with the same inputs and the resulting images looks as good as can be expected: `topup --verbose --imain=bud --datain=acq_params.txt --config=$FSLDIR/src/topup/flirtsch/b02b0.cnf --out=rs_topup --fout=topup_field

applytopup --imain=../bia6_00186_004_01.nii.gz --inindex=1 --method=jac --datain=func_params.txt --topup=rs_topup --out=run004 --verbose `

Screen Shot 2022-08-17 at 9 49 33 AM

If i look at my topup field side-by-side with the fmriprep version, it appears the fmriprep is being stretched and potentially flipped somewhere:

Screen Shot 2022-08-18 at 9 12 18 AM

As a sanity check, i went back to v20.#.# where AFNI could be used for pepolar correction and the results seem reasonable with no other variations to my fmriprep call, except the -max-topup-vols:

Screen Shot 2022-08-18 at 9 19 29 AM

Just looking for some guidance on where the fieldmap correction may be going wrong and how to get a clean run.

oesteban commented 2 years ago

Thanks much for reporting @cmpetty.

May I ask you to re-run this participant with the argument --debug fieldmaps added? Then, please double check the fieldmap report (the one you pasted under The actual SDC correction looks reasonable on the pepolar images:).

You're right it seems TOPUP is performing properly. Let's now see whether SDCFlows is doing an appropriate handling of the fieldmap.

oesteban commented 2 years ago

Sorry, I'm still digesting all the info in your post (thanks for the thoroughness, it's really useful!)

If i look at my topup field side-by-side with the fmriprep version, it appears the fmriprep is being stretched and potentially flipped somewhere:

I would not expect any meaningful change, but just to be sure: could you regenerate one of the figures at the same z-slice?

cmpetty commented 2 years ago

i ran v22 with the debug fieldmaps option, but i didn't see any debug info in the report

I uploaded the results and my own topup run here: https://github.com/cmpetty/fmriprep_debugging

A couple observations i noticed when debugging.

When i ran topup with the merged fieldmaps that fmriprep produced, my results looked identical ( poor ). fmriprep merged them as PA, AP. The functionals were acquired AP.

Screen Shot 2022-09-02 at 10 22 49 AM

When i re-ordered the pepolar images to be AP, PA ( like my original post ) my topup results looked better: Screen Shot 2022-09-02 at 3 14 09 PM

A side-by-side of my topup_field and fmriprep's "corrected field": Screen Shot 2022-09-02 at 10 19 54 AM

A side-by-side of my topup field and the results from topup, using the PA/AP ordered blips: Screen Shot 2022-09-02 at 3 31 58 PM

topup --verbose --imain=sub-3152_ses-intake_acq-orig_dir-PA_epi_merged_volreg.nii.gz --datain=acq_params_fmrip.txt --config=$FSLDIR/src/topup/flirtsch/b02b0.cnf --out=rs_topup_fmrip --fout=topup_field_fmrip

A side-by-side of the applied topup field AP/PA ( mine ) and PA/AP ( fmriprep ) on the mean image:

Screen Shot 2022-09-02 at 3 45 40 PM

the 2nd one is consistent with the fmriprep v21+ output and the first is what i would've expected.

julfou81 commented 2 years ago

Very interesting case and discussion. Thank you for all this reporting. I don't have any solution to propose but rather some thinking about all this. The different behavior of topup depending of the order of the AP/PA merge may have to do with the fact that the first volume of this merge will define the space in which the distorsion field will be calculated. Also, in such a strong distorsion case, it could be that the regular topup parameters (defined in the b02b0.cnf file) are not optimal and the topup calculation stopped too soon. I guess that if you look at the corrected AP images, they must be still quite different from the corrected PA images at the end of topup execution. What could be the solution for you? Using fmriprep, as it stands, as you noticed in such extreme case, 3dQwarp from AFNI may give better results (with standard parameters in both case). One solution could be to feed FMRIPREP v21+ with specific topup configuration file when needed but it is not possible as of now. You could try the Fieldmapless option but I am pessimistic on this, from personal experience it did not work well for cases with similar extreme distortions (also a mistake in the acquisition with of shimming box much too small).

oesteban commented 2 years ago

The different behavior of topup depending of the order of the AP/PA merge may have to do with the fact that the first volume of this merge will define the space in which the distorsion field will be calculated.

Indeed this was one thing we implemented in nipreps/sdcflows#278 (towards addressing #2830). Let's see if that was the culprit.

Also, in such a strong distorsion case, it could be that the regular topup parameters (defined in the b02b0.cnf file) are not optimal and the topup calculation stopped too soon. I guess that if you look at the corrected AP images, they must be still quite different from the corrected PA images at the end of topup execution.

SDCFlows does not use the default parameters, and defines its own (https://github.com/nipreps/sdcflows/blob/6061ed8a215d4be4b77c253f647343244113f2ee/sdcflows/data/flirtsch/b02b0.cnf). However, I do agree that this is a line of future investigation/optimization.

What could be the solution for you? Using fmriprep, as it stands, as you noticed in such extreme case, 3dQwarp from AFNI may give better results (with standard parameters in both case).

As mentioned in #2830, I strongly recommend against the 3dQwarp idea (unless a new version addressing the limitations is released).

One solution could be to feed FMRIPREP v21+ with specific topup configuration file when needed but it is not possible as of now.

This could be investigated, but I would first make sure that this is not just a problem of SDCFlows getting the PE polarities wrong at some point. We will cut a release candidate soon so that @cmpetty can test with the latest fixes.

You could try the Fieldmapless option but I am pessimistic on this, from personal experience it did not work well for cases with similar extreme distortions (also a mistake in the acquisition with of shimming box much too small).

That option should not be used unless you don't have any fieldmaps. We need to fix this issue :)

Thanks for the useful comments! Happy to get you more involved in SDCFlows, @julfou81 -- we really need contributors who know a lot about susceptibility distortions, topup, etc.

Ludoalevesque commented 1 year ago

Hi, We are also experimenting a similar issue with our fMRI data. We have some heavily distorted data for some subjects. Here is an example of the raw data for both phase directions (the fMRI was acquired PA) :

sub-10_bold_raw sub-10_rev_raw

After running fmriprep we had weird results for the SDC. Here is the bold ref and the preprocessed bold as outputted by fmriprep :

sub-10_fmriprep_orig_boldref sub-10_fmriprep_orig_preproc-bold

We also noticed that the fmap outputs were weird, as you can see in the fmap coefficient image :

sub-10_fmriprep_orig_fmap_coef

We were also told this might have something to do with skull striping, so we tried to run bet on the both the BOLD and the reverse images before fmriprep. Here is what it did. Fieldmap coefficient:

sub-10_fmriprep_bet-both_fmap_coef

This image seemed better than what we had previously, but the bold ref and the preprocessed bold were still showing some weird artifacts. Bold ref:

sub-10_fmriprep_bet-both_boldref

Preprocessed bold:

sub-10_fmriprep_bet-both_preproc-bold

Considering that functional inputs for fmriprep should not normally be processed data, we tried to run fmriprep on the raw bold, but in the fmap folder, we kept the brain extracted reverse epi, to see if the results were better.

sub-10_fmriprep_bet-rev_fmap_coef sub-10_fmriprep_bet-rev_boldref sub-10_fmriprep_bet-rev_preproc-bold

Surprisingly, the results look a bit better in the prefrontal cortex, but we don’t really understand why …

We also ran topup on our own and the results also seemed better:

sub-10_topup_coef sub-10_topup_corrected

*this one looks a bit different from the other fmriprep outputs, but keep in mind that this is in native space while fmriprep outputs are not.

We don’t really understand what is happening. We hope that brings more information to the issue and that someone can give us more insight on where this problem might be coming from. Thank you !

effigies commented 1 year ago

This seems to be the same as https://github.com/nipreps/fmriprep/issues/3013. The issue appears to be with the Jacobian correction, which we have been omitting. We need to resolve this soon, because these distortions are part of the correction process, and they need to be attenuated by applying a Jacobian weighting.

effigies commented 11 months ago

This should be resolved in 23.2.0a1.