How do I know what to put into my --acqp file?

The --acqp file is of the same format as the --datain parameter for topup, and quite possibly you can re-use the one you created for topup. It is used to inform eddy what direction the distortions are likely to go in. Each row consists of a vector (three values) that specify what the phase-encode (PE) axis is and also what direction along that axis that imply higher frequency. So for example the two lines

0 1 0 0.085
0 -1 0 0.085

have the vectors [0 1 0] and [0 -1 0] which both imply that PE is along the y-direction and which further imply that for the first row [0 1 0] a higher frequency is associated with a position higher up along the y-direction (i.e. positive blips) and the second [0 -1 0] implies that a lower frequency is associated with a position higher up (i.e. negative blips).

The fourth element in each row is the time (in seconds) between reading the center of the first echo and reading the center of the last echo. It is the "dwell time" multiplied by "number of PE steps - 1" and it is also the reciprocal of the PE bandwidth/pixel.

For a Siemens scanner you will get a "protocol PDF" that contains all the relevant information. Look for the tags "Phase enc. dir.", "Echo spacing" and "EPI factor". Here is an example

Phase enc. dir. A >> P
...
Echo spacing 0.8 [ms]
...
EPI factor 128

Where the first tag tells us that more posterior (lower down in the y-direction) locations are associated with higher frequency meaning that our vector should be [0 -1 0] and where the second two tags tells us that the fourth number should be 0.8*0.001*127=0.102 resulting in

0 -1 0 0.102

as the row for our --acqp file.

I am afraid that I don't have sufficient experience of the other scanner manufacturers to know how/where the relevant information is stored.

Does it matter what I put into my --acqp file?

Almost not at all, as long as the same file is used for both topup and eddy.

Let us for example say that the "true" file should be

0 1 0 0.05
0 -1 0 0.05

but that we happen to use one that looks like

0 -1 0 0.1
0 1 0 0.1

instead. Surely that must be disastrous? No, in fact the results would be indistinguishable from each other.

The sign-swap for the PE-directions would mean that topup would sign-swap the entire field so that areas with higher field will look like they have a lower field, and vice versa. The twice as large "readout time" means that topup will underestimate the field by a factor 2. Hence, the estimated field will be sign-swapped and underestimated, which on the surface sounds like a bad thing.

But now we enter that field into eddy with the --topup parameter and also the "fault" --acqp file. The sign-swap in the field and in the --acqp file will now cancel out by being multiplied together, so the displacements will still go in the correct direction. Similarily, the field is underestimated by a factor 2, but that is counteracted by a twice as large "sensitivity" to off-resonance indicated by the too large readout time. So the magnitude of the displacements is also correct.

There are some special cases where it matters to get the --acqp file right, but unless you know exactly what you are doing it is generally best to avoid those cases. They would be

If your data is not one of those special cases, and you don't have a specific interest in the off-resonance field per se, you can make life very easy for you. Run a movie (for example using fslview or FSLeyes) of the 4D file you plan to enter into topup. Does the brain jump up and down? If so, use

0 1 0 0.05
0 -1 0 0.05

Does the brain bounce from side to side? If so, use

1 0 0 0.05
-1 0 0 0.05

Is the brain essentially still? If so, there has probably been a mistake in the acquisition and both images have been acquired with the same phase-encoding. In that case you cannot use topup.

What if I still really want to have an "accurate" --acqp file?

If you are one of the special cases, or you just really would like to see what the "true" off-resonance field is, here is an --acqp for dummies:

What you see in FSLVIEW

P_2_A.jpg

A_2_P.jpg

R_2_L.jpg

L_2_R.jpg

In Siemens protocol

Phase enc. dir. P >> A
Echo spacing 0.75 [ms]
EPI factor 128

Phase enc. dir. A >> P
Echo spacing 0.75 [ms]
EPI factor 128

Phase enc. dir. R >> L
Echo spacing 0.96 [ms]
EPI factor 128

Phase enc. dir. L >> R
Echo spacing 0.96 [ms]
EPI factor 128

In --acqp file

0 1 0 0.095

0 -1 0 0.095

1 0 0 0.122

-1 0 0 0.122

Do I ever need more than two rows (for eddy) in my --acqp file?

I would say very rarely, if ever.

The rationale for it was to use the movement estimates from the topup output as starting estimates for eddy's movement estimation. Changes to the details of how eddy performs the movement estimation has made that less important. As of the FSL 5.0.9. eddy patch, eddy will start by rigid-body registering all the b=0 and then linearly interpolate the movement parameters from those to the diffusion weighted images. If the b=0 volumes are interspersed among the diffusion images this works really well even in the presence of large subject movements. If additional help is still needed one can typically achieve that by adding --niter=8 and --fwhm=10,6,4,2,0,0,0,0 to the eddy command line.

Having said that, there is nothing wrong with having more than two rows in the --acqp file, and it can certainly do no harm.

Will eddy rotate my bvecs for me?

Yes, it will as of the FSL 5.0.9 patch and later versions. Note that for the rotation to work the bvecs fed into eddy must be correctly oriented for use by the FSL tools.

What would a good eddy command look like for data with lots of movement?

There are several things we need to take into consideration when there is lots of movement:

One example of an eddy command that does all of that would be

eddy_cuda --imain=my_data --acqp=acqparams.txt --index=index.txt --mask=my_brain_mask --bvals=bvals --bvecs=bvecs --topup=my_topup --out=my_hifi_data --niter=8 --fwhm=10,6,4,2,0,0,0,0 --repol --ol_type=both --mporder=8 --s2v_niter=8 --slspec=my_slspec.txt

In this example --niter=8 and --fwhm=10,6,4,2,0,0,0,0 instructs eddy to run 8 initial iterations with decreasing amounts of smoothing. During these iterations eddy estimates eddy current-induced fields, inter-volume movement parameters and detects and replaces outliers. The outlier detection is requested with the --repol flag. The --ol_type=both parameter instructs eddy to detect both groups of outliers (in a multi-band acquisition, which we assume this data set is, all slices in an MB-group are affected in the same way by subject movement) and individual outlier slices (pulsatile movement of the brain can cause one slice in an MB-group to be an outlier, but not affect the others).

The --mporder=8 and --s2v_niter=8 parameters tells eddy to run an additional 8 iterations. For the additional iterations the movement model is replaced by a slice-to-vol method with 9 (8+1) degrees of freedom per movement parameter and volume. During these additional iterations the estimates of the eddy current-induced fields and the outliers are also further refined. The parameter --slspec=my_slspec.txt specifies a text file with information of acquisition order and multi-band structure. You can find more information about the --slspec file here and here.

Note also that I have specifically written eddy_cuda. That is because the slice-to-vol movement model is only implemented for the CUDA version.

How should my --slspec file look?

The --slspec file is necessary when using the slice-to-vol movement correction introduced in version 5.0.11. It is used to inform eddy of the temporal order in which the slices of a volume was acquired, and also to specify how multiple slices were grouped in a multi-band (MB) acquisition. In certain cases it is quite easy to work this out, for example in the case of a single band acquisition or in the case of an MB acquisition with an odd number of excitations/groups when using the UMinn MB sequence on a Siemens scanner. But even in the case of a single band acquisition there are some "gotchas" that one need to be aware of.

Hence, when in doubt it is always safest to use the information in the DICOM headers to create the --slspec file. This can be done either directly from the DICOM files, or via a .json file created by your DICOM->niftii conversion software. If you for example have the following file obtained from a conversion with dcm2niix

Sample JSON file from dcm2niix

mb_diff_1007_V1_10.json

you can then create your --slspec file using the following Matlab code

fp = fopen('name_of_my_json_file.json','r');
fcont = fread(fp);
fclose(fp);
cfcont = char(fcont');
i1 = strfind(cfcont,'SliceTiming');
i2 = strfind(cfcont(i1:end),'[');
i3 = strfind(cfcont((i1+i2):end),']');
cslicetimes = cfcont((i1+i2+1):(i1+i2+i3-2));
slicetimes = textscan(cslicetimes,'%f','Delimiter',',');
[sortedslicetimes,sindx] = sort(slicetimes{1});
mb = length(sortedslicetimes) / length(find(diff(sortedslicetimes)==0));
slspec = reshape(sindx,[mb length(sindx)/mb])'-1;
dlmwrite('my_slspec.txt',slspec,'delimiter',' ','precision','%3d');

No doubt there are more elegant ways of doing it, but you would have to write those yourselves.

I look at my eddy-corrected data, and it looks like there is a mis-registration between the b=0 and the dwi volumes.

First of all it should be said that it is not easy to tell by just eye-balling the data. The high signal from CSF in the b=0 data, along with a quite different distribution of CSF on different sides of the brain tends to "trick the eye". What I recommend doing when checking the alignment is to focus on points of high curvature (typically at the "bottom" of sulci) that can be identified in both b=0 and dwi images.

Having said that, it is indeed possible to have such a mismatch. Inside eddy the dwi:s and the b=0 volumes are registered separately. In the initial version of eddy the first b=0 and the first dwi volume were assumed to be in register (by virtue of having been acquired closely in time) and that was that. Later versions have attempted to be a little more clever about it. As of eddy version 5.0.11 we hope to have sorted this out. There is a long description of that here and also how you can affect how it is done through the --dont_sep_offs_move and --dont_peas flags.

What does it mean if I have more than two rows in my --acqp file?

First of all, all the information below is correct. But as of the release of the FSL 5.0.9. eddy patch (here) it is also a little less relevant. The reason for that is that as of that version eddy will start by rigid-body registering all the b=0 and then linearly interpolate the movement parameters from those to the diffusion weighted images. If the b=0 volumes are interspersed among the diffusion images this works really well even in the presence of large subject movements. Hence, it is now quite rare that eddy needs the kind of help that the multiple --acqp rows offer.

Let us say you have data that has only been collected in two different ways (e.g. R->L and L->R) but your --datain file from topup contains four rows. Why do I need four rows? Should not two rows be sufficient to describe the two ways that my data was collected in?

The reason that more than two rows can be useful is that it is a way of feeding information from topup to eddy. Let us say we have collected two b=0 and two dwis each for R->L and L->R (which is of course silly, but used here to save space)

eight_original_images_b0_dwi_b0_dwi_b0_dwi_b0_dwi.jpg

Let us further say you used fslroi and fslmerge to create a file my_b0.nii.gz with only the b=0 scans

four_original_images_b0_b0_b0_b0.jpg

In order to run topup on this we create a --datain file (that we call parameters.txt) which consists of

-1 0 0 0.051
-1 0 0 0.051
1 0 0 0.051
1 0 0 0.051

since topup expects one line per volume. The output from topup consists of an estimate of the susceptibility induced field that caused the distortions and also an estimate of the position of each of the b=0 scans relative the first one. It is necessary for topup to simultaneously estimate the field and the subject movement as the latter may otherwise affect the estimation of the field. If the --out parameter in the topup call was set to for example --out=my_topup_output there will be a text-file named my_topup_output_movpar.txt which will contain something like

0 0 0 0 0 0
0.725 -0.023 -0.066 0.002 0.000 -0.002
0 -0.111 -0.327 0.002 0.013 -0.004
-0.698 -0.120 -0.426 0.002 0.014 -0.004

which are the rigid-body parameters for each of the scans relative the first.

If we now consider using a two line --acqp when running eddy on the full data

eight_original_images_b0_dwi_b0_dwi_b0_dwi_b0_dwi.jpg

we would have

--acqp
-1 0 0 0.051
1 0 0 0.051

and


--index
1 1 1 1 2 2 2 2

This would be sufficient to inform eddy about the acquisition parameters but we have lost the information about the relative positions of the b=0 scans. This information is potentially useful as a starting estimate for eddy when it calculates the relative position of the diffusion weighted volumes. The thinking here is that if you for example acquire a b=0 volume every ten volumes then the position of the first diffusion weighted volume (coming immediately after the first b=0 volume) is likely to be similar to first b=0 volume. Likewise the position of the tenth diffusion weighted volume (coming immediately after the second b=0 volume) is likely to be more similar to that of the second b=0 volume. So, in order to pass this information (for our silly eight scan example) to eddy one can instead use the --acqp - --index pair

--acqp
-1 0 0 0.051
-1 0 0 0.051
1 0 0 0.051
1 0 0 0.051

and


--index
1 1 2 2 3 3 4 4

which carries the exact same information pertaining to the acquisition parameters as the previous (two line) pair, but which now also carries information about the relative positions by virtue of the --index file also indexing the my_topup_output_movpar.txt seen above.

 

eddy/Faq (last edited 13:08:41 31-08-2017 by JesperAndersson)