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

[cmpetty]

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:

The actual SDC correction looks reasonable on the pepolar images:

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

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
`

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:

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:

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

[oesteban]

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]

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]

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.

When i re-ordered the pepolar images to be AP, PA ( like my original post ) my topup results looked better:

A side-by-side of my topup_field and fmriprep's "corrected field":

A side-by-side of my topup field and the results from topup, using the PA/AP ordered blips:

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:

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

[julfou81]

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]

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]

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) :

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

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

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:

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:

Preprocessed 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.

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:

*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]

This seems to be the same as #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]

This should be resolved in 23.2.0a1.