Eddy_quad works with slspec, not header info

m.bart · May 6, 2021, 1:48pm

Hi Rob,

incorporated eddy_quad is perfect! Works fine when feeding dwifslpreproc with slspec file.

I am not sure how it was implemented and how well can it extract info from the header of the file. As I had dicom i directly converted it info .mif with slicetiming info within .mih part of .mif, I tried to not provide dwifslpreproc with slspec (I though It will find the information about Slice Time on its own). Dwifslpreproc run OK with creating quad folder and all of its usual content. But folder with qc.json and qc.pdf was missing. It works fine on the same data when I fed algorithm with SliceTime.txt file in slspec. Is it how it is supposed to work and text file should be always included or is it a mistake? Thanks for the answer.

Michaela.

rsmith · May 7, 2021, 6:07am

Hi Michaela,

If -eddyqc_text or -eddyqc_all are specified, and EddyQC is found on your system, it will attempt to run it; if just EddyQC fails, then dwifslpreproc will not abort at that point, it’ll still run to completion, but it should be issuing you with a warning at the command-line to say that EddyQC failed and therefore your “EddyQC output folder” will only contain those files output by eddy itself.

As per the documentation, dwifslpreproc will use slice timing information from the image header if it is available, such that you do not need to specify --slspec every time. The fact that EddyQC works when you provide that file manually, but not when reading the slice timing information from the header, means that the resultant slice group data being provided to EddyQC is different between those two cases. If you run dwifslpreproc with the -nocleanup option, you’ll be able to recover the slspec file that was generated internally by dwifslpreproc, and compare it against your manual version. There will also be a text file providing the error message that was generated by EddyQC.

My leading hypothesis is that your DICOMs report slightly different times for those slices within the same multi-band group, such that my gorgeous one-liner that converts slice times to slice groups does not yield the correct result. I have observed something similar-but-not-quite-the-same in my own data, where two DWI series acquired with everything equivalent except for reversal of the phase encoding direction could have fractionally different slice acquisition times reported, and hence attempting to concatenate those data would cause problems without explicit handling. There’s some chance that your problem may be cause by a similar effect, but due to minute variations in reported slice timings between members of the same group. If that’s the case, I’d need to modify dwifslpreproc to better handle such, but I can also give some advice w.r.t. MR protocol parameters that should stop it from appearing in your data to begin with.

Also, just in case, check if you have both multi-band and an odd number of slices, in which case you’re encountering this problem. Though I’d have expected in that case that manually providing the slspec file would encounter the same issue, since the relevant handling is upstream of the image padding…

Rob

m.bart · May 17, 2021, 10:11am

Hi rob,
I tested this issue out again and find out following: as you mentioned, all approaches run into same warning resulting with no Eddy QC .pdf and .json report produced (sorry for mystification). I am sure it is not a problem of FSL installation since on a 3 different dataset I did not faced this problem.

my dataset:

b_value (No of dirs): 0 (11), 750 (17), 1500 (26), 3000 (35), 5000 (43)
full rPE: AP_PA with identical dw_schemes (checked all 12 decimal places)
Dimension 120x120x75x132

I am giving here a prtsc from Thijs Dhollander’s script for generating image of b vectors distribution:

When I run dirstat I’ve got following warning for each shell with value varying from 0.14 to 0.26

I applied on a data same pipeline as usual, ( mrconvert dcm to mif; no denoising since it produces problems when eddy called too - we discussed it in a separate thread; mrdeggibs, dwifslpreproc). For the last step I have run following command:

dwifslpreproc ${i}_dwi_den_unr_AP_PA.mif ${i}_dwi_den_unr_preproc.mif -pe_dir AP -rpe_all -nocleanup -eddy_options "  --slm=linear --repol --mporder=16 --cnr_maps " -eddyqc_all  ${Dir}/${i}/Mrtrix3_proc/QC

During the dwifslpreproc I’ve got warning about not providing report but I am not able to uncover why:

Versions tested (slspec files copied from -nocleanup):

mrconvert .dcm to .mif, no slspec file provided

0 25 50
2 27 52
4 29 54
6 31 56
8 33 58
10 35 60
12 37 62
14 39 64
16 41 66
18 43 68
20 45 70
22 47 72
24 49 74 75
1 26 51
3 28 53
5 30 55
7 32 57
9 34 59
11 36 61
13 38 63
15 40 65
17 42 67
19 44 69
21 46 71
23 48 73

mrconvert .dcm to .mif, slspec file provided

0 41 57
25 66 9
50 18 34
2 43 59
27 68 11
52 20 36
4 45 61
29 70 13
54 22 38
6 47 63
31 72 15
56 24 40
8 49 65
33 74 17 75
58 1 42
10 26 67
35 51 19
60 3 44
12 28 69
37 53 21
62 5 46
14 30 71
39 55 23
64 7 48
16 32 73

provided by me was:

mrconvert .dcm to .nii to .mif (to have no info about slicetiming in a header of .mif), slspec file provided

0 41 57
25 66 9
50 18 34
2 43 59
27 68 11
52 20 36
4 45 61
29 70 13
54 22 38
6 47 63
31 72 15
56 24 40
8 49 65
33 74 17 75
58 1 42
10 26 67
35 51 19
60 3 44
12 28 69
37 53 21
62 5 46
14 30 71
39 55 23
64 7 48
16 32 73

Provided by me same as in a 2nd case.

Slspec files provided by me were generated with my own script. Original slicetiming from dicom header was as follows:

AP:
0.0000,1.6700,0.1300,1.7975,0.2575,1.9275,0.3850,2.0550,0.5150,2.1825,0.6425,2.3125,0.7725,2.4400,0.9000,2.5675,1.0275,2.6975,1.1575,2.8250,1.2850,2.9550,1.4125,3.0825,1.5425,0.0000,1.6700,0.1300,1.7975,0.2575,1.9275,0.3850,2.0550,0.5150,2.1825,0.6425,2.3125,0.7725,2.4400,0.9000,2.5675,1.0275,2.6975,1.1575,2.8250,1.2850,2.9550,1.4125,3.0825,1.5425,0.0000,1.6700,0.1300,1.7975,0.2575,1.9275,0.3850,2.0550,0.5150,2.1825,0.6425,2.3125,0.7725,2.4400,0.9000,2.5675,1.0275,2.6975,1.1575,2.8250,1.2850,2.9550,1.4125,3.0825,1.5425

PA:
0.0000,1.6700,0.1300,1.7975,0.2575,1.9275,0.3850,2.0550,0.5150,2.1825,0.6425,2.3125,0.7700,2.4400,0.9000,2.5675,1.0275,2.6975,1.1550,2.8250,1.2850,2.9525,1.4125,3.0825,1.5400,0.0000,1.6700,0.1300,1.7975,0.2575,1.9275,0.3850,2.0550,0.5150,2.1825,0.6425,2.3125,0.7700,2.4400,0.9000,2.5675,1.0275,2.6975,1.1550,2.8250,1.2850,2.9525,1.4125,3.0825,1.5400,0.0000,1.6700,0.1300,1.7975,0.2575,1.9275,0.3850,2.0550,0.5150,2.1825,0.6425,2.3125,0.7700,2.4400,0.9000,2.5675,1.0275,2.6975,1.1550,2.8250,1.2850,2.9525,1.4125,3.0825,1.5400

It took a while to notice but as you mentioned, there are small differences between AP and PA acquisition.

During writing this answer, my eyes dropped on an odd number of slices I have here and never had before . I run into this answer from Jesper Andersson on JISCMail. He suggests not to delete any slices (as I used to before - because of multiband), but to use a different config file. So I have specified it within dwifslpreproc as `-topup_options " --config=b02b0_1.cnf " Although according to FSL web (topup/TopupUsersGuide - FslWiki), I don’t need to redirect to FSLDIR directory, it searches within automatically, It didn’t (it is just a minor advise for future readers). Answer from terminal when command run:

It looks like command was called with 2 config file arguments? From original 10 hours, time for Topup step increased to 91.6 hours only to get the same warning and not producing files again.

I also have one minor question. I couldn’t help but notice that the slspec file from -nocleanup (provided by dwifslpreproc) when using only .mif with SliceTiming info directly from .dcm looks different then mine from a script I made. My matlab code (MBF stands for multi band factor, here it is 3). Data were measured as interleaved on Siemens Prisma 3T:

[~, slice_order] = sort(info.dicom_header.SliceTiming);
slice_order_reshape=reshape(slice_order,length(slice_order)/MBF,MBF);
slice_order_reshape=slice_order_reshape-1;

Thanks a lot for your time and any answers!

Stay happy and healthy,
Michaela.

rsmith · June 25, 2021, 7:55am

Hi Michaela,

It seems you have stumbled inadvertently into my own little slice of hell.

Having neither approach work makes a lot more sense, especially in light of showing your slspec files. What you’ve uncovered is the issue described here, for which I currently have a proposed solution sitting here that I unfortunately don’t currently have the capability to easily test myself. If you’ve the time / motivation to compile from source and try it out on your data, please do.

In a much earlier version of what is now dwifslpreproc, topup would crash if provided with an image with an odd number of voxels in any spatial axis, and FSL only came with a single configuration file. So one of the “capabilities” that this script would provide over and above being a dumb-and-blind execution of topup & eddy was to pad the image data, run these commands, then remove the padding. Then slice-to-volume registration was added, which requires the slspec file; so if it was specifically the number of slices that was odd in the input data, I had to add that padded slice to the corresponding slice group (note that the slspec padding needs to echo the image padding regardless of where the slspec data came from originally, which is why the second table you provided differs from the third. This worked, until EddyQC was released, which doesn’t like slspec files that aren’t square matrices (which I argued should be entirely permissible, but lost). If you were to dig through the script scratch directory to find the EddyQC output text file, it should contain some cryptic Python error about a non-square matrix. So what the proposed solution does is remove the padding that the dwifslpreproc script added to both the image and slspec data, including everything that was generated by eddy, prior to invoking EddyQC.

During writing this answer, my eyes dropped on an odd number of slices I have here and never had before … From original 10 hours, time for Topup step increased to 91.6 hours

While I did see that alternative topup config files with different downsampling ratios were added in I think FSL 6.0.0, I knew that without image downsampling the execution times would be even longer than those that users are already complaining about, so didn’t want to rely exclusively on that approach. This pad - unpad approach is just a different way of dealing with the fact that when topup downsamples, it strictly demands that the downsample factor be a divisor of the voxel count in all axes.

It looks like command was called with 2 config file arguments?

Correct: dwifslpreproc will always insert the --config= entry into the function call, so if you also specify it yourself in -topup_options, it will appear twice. Preventing this would require instructing dwifslpreproc to search for the presence of “--config=” in the content of the input to -topup_options (and/or via a new command-line option -topup_config), explicitly copy that config file into the scratch directory, and instruct topup to use the config file at that location instead of the default one, which is something I’ve not yet done. Exactly what topup will do when it sees that option appear twice, I don’t know: MRtrix3 commands will error instantly if a command-line option that can only meaningfully be specified once appears multiple times, but topup is outside of my control.

It took a while to notice but as you mentioned, there are small differences between AP and PA acquisition.

Yep, exactly the same as I’d observed myself; but isn’t the source of the issue as it turns out, it was just my naive hypothesis.

I couldn’t help but notice that the slspec file from -nocleanup (provided by dwifslpreproc) when using only .mif with SliceTiming info directly from .dcm looks different then mine from a script I made.

I’m always hesitant about inferring such data from “expectations” about the acquisition rather than sidecar data. But given that the first column in the first three rows corresponds to the first slice group, I’d suggest you’re missing a transpose operation somewhere.

Cheers
Rob