Qsiprep - Distortion correction not applied to eddy_corrected image?

Summary of what happened:

Hello!

I am running qsiprep preprocessing on dwi data acquired with reversed phase-encoding directions (AP and PA), each with 93 diffusion directions. It completes successfully. My topup results look good (are corrected for susceptibility distortions), but my eddy results don’t appear to be corrected for these distortions, and likely neither the eddy current distortions. See screenshots below.

I have tried running qsiprep both with and without the --distortion-group-merge flag (set to average) and get similar results. What am I missing?

Command used (and if a helper script was used, a link to the helper script or the command generated):

docker run -ti --rm --gpus all \
	-v /usr/local/freesurfer/7.4.1/license.txt:/opt/freesurfer/license.txt:ro \
	-v $bids_dir:/data:ro \
	-v $eddy_config_file:/sngl/eddy/eddy_config.json:ro \
	-v $out_dir:/out \
	-v $out_dir:/scratch \
	pennbbl/qsiprep:0.20.0 /data /out participant \
	--eddy-config /sngl/eddy/eddy_config.json \
	--output-resolution $res --pepolar-method TOPUP \
	--distortion-group-merge average \
	--bids-database-dir $bids_dir \
	-w /scratch 

Version:

qsiprep 0.20.0

Environment (Docker, Singularity / Apptainer, custom installation):

Docker

Data formatted according to a validatable standard? Please provide the output of the validator:

PASTE VALIDATOR OUTPUT HERE

Relevant log outputs (up to 20 lines):

PASTE LOG OUTPUT HERE

Screenshots / relevant information:

eddy_corrected image from the qsiprep_wf/single_subject_0001_wf/dwi_preproc_wf/hmc_sdc_wf/eddy directory:

topup_imain_corrected image from the qsiprep_wf/single_subject_0001_wf/dwi_preproc_wf/hmc_sdc_wf/topup directory:


Hi @jhau,

It would probably be better to look at the quality of the final preprocessed images, including screenshots in the html.

Best,
Steven

Hi @Steven,

From the html, here is a screenshot of the topup results:

and a screenshot of the b0 ref image before correction, which already looks odd - the after image is identical just cropped to the blue outline:

Thanks for your help.

Best,
Jan

It looks to me like there might be a problem with the PhaseEncodingDirection in the sidecars. Are there separate AP/PA images in your fmap directory? It’s possible that topup is working, but the flipped PE Direction for the dwi images is causing the wrong correction to be applied during eddy.

Hi @mattcieslak,
I have my dir-AP and dir-PA images, labelled correspondingly, in the dwi directory. No fmap directory as I don’t think that’s needed in this case, correct? The json files show the PhaseEncodingDirection field to be j- for the AP and j for the PA images, which I believe are correct.

As an additional check, I tried swapping the AP and PA labels (i.e., I renamed the AP nii/bval/bvec/json files to PA and PA nii/bval/bvec/json files to AP) and reran qsiprep. I get the same results as before - the topup results look good, but the b0 ref is not uncorrected for distortions. Based on this, I don’t think the results are due to my PE Directions being flipped. Any other ideas?

Following up on this. Has anyone encountered this type of issue with eddy before?

My topup results look great but my eddy_corrected image looks bad (oddly stretched and shows some ghosting which isn’t present in the raw data and actually looks worse than the raw data).

I can’t figure out what is amiss in my eddy step. I’ve checked the eddy_acqp.txt and eddy_index.txt and those correspond with my input image (merge__merged) and the brain mask looks fine. I also ran eddy swapping the PE of the rows in my eddy_acqp.txt but it gives even worse results. I know it’s doing some sort of susceptibility distortion correction but why am I not getting results in my dwi image that look like the topup corrected b0? Maybe something is off in the estimation of the eddy current distortions?

Hi @jhau, could you please try with the most recent version (0.21.4)? We fixed two topup-related bugs that could sometimes result in something like this.

@mattcieslak, I did run it with the latest version QSIPrep 0.21.5.dev0+g36b93fe.d20240507 and still get the same results unfortunately.

Are you using a new working directory each time?

Yes, I’m using a new output and working directory each time.

I ran eddy using the --topup flag specifying the field coefficient and movement parameters and that produced good results (matches topup_imain_corrected). The issue I am having is with the --field and --field_mat input. FSL notes some caveats when using --field and recommends using --topup. I’ve triple-checked my inputs to eddy (acqp/index settings and concatenated dwi file) and they should be correct so I don’t know what is causing the --field and --field_mat to give poor results. Could we set eddy to use the --topup flag instead of --field when using the topup method to avoid this?

Command:
eddy_cuda10.2 --ol_nstd=4 --ol_type=both --cnr_maps --estimate_move_by_susceptibility --imain=qsiprep_wf/single_subject_0426_wf/dwi_preproc_wf/pre_hmc_wf/rpe_concat/merge__merged.nii.gz --index=qsiprep_wf/single_subject_0426_wf/dwi_preproc_wf/hmc_sdc_wf/gather_inputs/eddy_index.txt --mask=qsiprep_wf/single_subject_0426_wf/dwi_preproc_wf/hmc_sdc_wf/transform_mask_to_eddy/topup_imain_corrected_avg_trans_mask_trans_flirt.nii.gz --acqp=qsiprep_wf/single_subject_0426_wf/dwi_preproc_wf/hmc_sdc_wf/gather_inputs/eddy_acqp.txt --bvals=qsiprep_wf/single_subject_0426_wf/dwi_preproc_wf/pre_hmc_wf/rpe_concat/merge__merged.bval --bvecs=qsiprep_wf/single_subject_0426_wf/dwi_preproc_wf/pre_hmc_wf/rpe_concat/merge__merged.bvec --out=qsiprep_wf/single_subject_0426_wf/dwi_preproc_wf/hmc_sdc_wf/eddy/eddy_corrected_topup_base --repol --residuals --json=qsiprep_wf/single_subject_0426_wf/dwi_preproc_wf/pre_hmc_wf/merge_plus/merge_dwis/merged_metadata.json --topup=qsiprep_wf/single_subject_0426_wf/dwi_preproc_wf/hmc_sdc_wf/topup/topup_imain_base

We don’t use the --topup flag because we select the “best b=0” images from each phase encoding direction, which may not be the first b=0 image in the dwi. Instead we run topup on the best b=0s and align the best b=0 (with the same PED) to the first image from the dwi scan sent to eddy.

Is this issue happening with all your data or just a few subjects?

Also, what are the orientations of your original fmap images and dwis?

I ran the same participant using qsiprep 0.16.1 and the sdc results look good:

I don’t think it’s specific to this subject but I will try on a few more to see if it’s reproducible.

I have a pair of dwi images with reversed phase-encoding (one is AP and the other is PA).

Nevermind, I spoke too soon. It is the same problem in version 0.16.1. The sdc (topup) results are good but for some reason the application of the fieldmap (estimated in topup) to eddy is bad.

I tested 3 subjects in qsiprep 0.20.0 and their preprocessed dwi images all have the same problem.

I know why, and this makes perfect sense. The TotalReadoutTime is the culprit*. When using the --field flag in eddy and the acqp file is different for topup and eddy, which in the qsiprep wf it is in my case for all 3 subjects that I tested (a result of the best b0 selection?), it’s very important that the metadata be accurate for correct application of the fieldmap. This isn’t an issue when the same acqp file is used in both topup and eddy.

Why is the acqp file different for eddy than topup in the qsiprep wf? If the best b0 is selected from each dwi image (AP and PA), that shouldn’t affect the order of AP and PA in topup/eddy, no?

*This data was acquired on a GE using the mux epi research sequence from Stanford and requires proprietary recon processing from the P-files. dcm2niix fails on this data but can still generate the metadata, but not reliably for all parameters. For instance, SliceTiming doesn’t account for the multi-band factor and I’m fairly certain now the TotalReadoutTime is off.

I ran eddy using the --topup flag (instead of --field) on my 3 subjects and all produced much better SDC results, if slightly blurry. This is due to the best b0s not being the first b0s in my dwi images as you pointed out and is why you don’t use the --topup flag.

P.S. If anyone has tips on calculating TotalReadoutTime for the GE mux epi sequence that would be greatly appreciated.

@jhau I’m glad you figured this out!! This definitely deserves a section in the documentation, as GE has been difficult for a lot of people

@mattcieslak, I am able to successfully correct for image distortions when using eddy with the --field flag in FSL, but the same subject fails in qsiprep. Any ideas?

Correction in my previous post: My acqp inputs are the same for both topup and eddy, just that they are ordered differently and in the eddy_index file it shows 2 2 2 … 1 1 1.

From the eddy documentation:
“If one uses the same --acqp file for both topup and eddy it doesn’t matter if one get the total acquisition time or the PE-polarity wrong. The errors in topup and eddy will cancel out and the end results will still be correct.”

I also ran qsiprep on one extracted AP_b0 image (using it as a fmap for my PA_dwi image) and get the same (bad) eddy results.

@mattcieslak, I believe I am having exactly this issue, which has a fix: FIX: Revise the TOPUP workflow and related utilities by oesteban · Pull Request #278 · nipreps/sdcflows · GitHub. My dwi data are oriented in LAS and appear to be one of these cases.

It looks like the qsiprep team are in the process of updating NiPreps. Looking forward to testing out the next release.