Unwarping BOLD

Part 6: BOLD Unwarping and Registration

Finishing basic preprocessing of BOLD data is complex, as the process of correcting distortions involves the highres anatomical scans, fieldmap scans, and the EPI scans. First, we will prepare the fieldmaps for use with distortion correction.

Preparing Fieldmaps

A fieldmap scan involves collecting two images with a slight delay between them. This allows us to estimate how much inhomogeneity there is in the local magnetic field at each voxel, which in turn provides an estimate of how much an EPI scan would have been distorted at different locations. This information is contained in the phase image from the fieldmap scan. The scan also provides a magnitude image, which has enough anatomical detail to register the fieldmap images to the highres anatomical.

The phase image provides an estimate of local changes in the fieldmap, but in order to obtain a global estimate to use for unwarping, the phase image must be "unwrapped" to obtain a radians/second image. FSL provides a tool for this operation. The fieldmap is noisy at the edges of the brain, so the maintainers of FSL recommend using a tight mask that excludes these noisy portions. To ensure we have a tight mask that excludes voxels outside the brain, we will use the mask we determined from the FreeSurfer output. This requires first registering the fieldmap magnitude image with the highres scan. Because distortion may vary over time as the participant's head moves, it is generally best to align to the highres scan taken closest in time to the fieldmap. In this case, we will align the day 1 fieldmap to the day 1 highres scan, and the day 2 fieldmap to the day 2 highres scan.

For example, to display the commands used to prepare fieldmap 1 for bender_03:

prep_fieldmap.py bender_03 -a 1 -f 1 --dry-run

The -a option indicates which anatomical scan to use; -a 1 indicates that we will align the fieldmap to orig_brain1.nii.gz. The -f option selects a fieldmap scan to use, in this case fieldmap_mag1.nii.gz and fieldmap_phase1.nii.gz. To submit a job to prepare both fieldmaps for both subjects, we can call slaunch twice, once for each fieldmap. In each case, we indicate which fieldmap to use and which anatomical image to align it to:

slaunch "prep_fieldmap.py {} -a 1 -f 1" bender_03:bender_04 -N 1 -n 4 -a 6 -r 00:05:00
slaunch "prep_fieldmap.py {} -a 2 -f 2" bender_03:bender_04 -N 1 -n 4 -a 6 -r 00:05:00

This will run the commands:

prep_fieldmap.py bender_03 -a 1 -f 1
prep_fieldmap.py bender_04 -a 1 -f 1
prep_fieldmap.py bender_03 -a 2 -f 2
prep_fieldmap.py bender_04 -a 2 -f 2

You should now have images in each subject's fieldmap directory named fieldmap_rads_brain1.nii.gz and fieldmap_rads_brain2.nii.gz.

Calculate Unwarping of EPI scans

This is a complicated step. We need to unwarp the EPI scan. The fieldmap provides an estimate of how much distortion there would have been at each voxel. But we need to align the EPI scan to the fieldmap in order to know what parts of the image were distorted and in which direction. But the alignment of the EPI scan to the fieldmap will be affected by the warping of the EPI. Therefore, alignment and unwarping must happen together. This process is handled by the FSL tool epi_reg.

In order to register the EPI to the fieldmap, the highres anatomical scan is used a common target for both the fieldmap and EPI. The EPI is aligned to the highres scan using a process called boundary-based registration. The white matter mask we created earlier is used to align the EPI scan to the anatomical scan, by maximizing the change in intensity of the EPI at points crossing the boundary between white matter and gray matter in the anatomical scan. This usually works quite well, but can cause some errors where the alignment "slips" and aligns the wrong intensity boundary in the EPI to the white matter boundary in the anatomical scan. Therefore, we will use epi_reg_ants, which starts with the registration estimate from epi_reg and refines it using a different, intensity-based registration method implemented in ANTS. This process is slow, since registration and unwarping must be done through multiple iterations, but accurate.

There are two things you need to know about your EPI scans before you can unwarp them.

• Effective echo spacing - this is the echo spacing in seconds (here, 0.00077), divided by the GRAPPA factor (generally 2 for our scans).
• Phase encoding direction - this is generally "y-" for our scans.

For a given functional scan, you must select an anatomical scan and a fieldmap scan to use. As mentioned before, if you have multiple days, select the scans taken on the same day as the functional you are processing.

For example, to see the commands that will be run for bender_03, prex_1 (this was taken on the first day):

epi_reg_run.py bender_03 prex_1 0.000385 -a 1 -f 1 -p y- -k --dry-run

The -k option indicates that we want to keep the intermediate files. If you're processing a lot of subjects, you will probably want to omit the -k option, so that the large intermediate files will be deleted to save space.

First, we'll process the prex_1 scans for both subjects:

slaunch "epi_reg_run.py {} prex_1 0.000385 -a 1 -f 1 -p y- -k" bender_03:bender_04 -N 1 -n 2 -a 12 -p development -r 00:40:00

Then the study_1 scans:

slaunch "epi_reg_run.py {} study_1 0.000385 -a 2 -f 2 -p y- -k" bender_03:bender_04 -N 1 -n 2 -a 12 -p development -r 00:40:00

Registration and Final Transformation

Next, we need to register the different functional runs to one another. We'll start from the unwarped average images that were calculated in the previous step, so that differences in warping between different scans have less of an effect on the co-registration.

You must first decide on which run will be the reference to provide a common functional scan space for each participant. Here, we'll use the study_1 run. We'll call reg_unwarp_bold_run.py, which processes one run.

Unlike some packages, this script uses nonlinear registration through ANTS. The advantage is that there may be distortion due to gradient nonlinearities that is different between different functional scans, as well as residual warping from local magnetic susceptibility gradients that was not perfectly corrected by the unwarping procedure; these nonlinear differences between scans can be corrected by nonlinear registration. There is a potential drawback to nonlinear registration: different functional scans may have dropout at different locations in the brain, and this may affect the registration. However, in tests there didn't seem to be large effects of differential dropout on the registration. This may be because ANTS registers images based on local cross-correlation, which should be relatively robust to intensity gradients. See below for discussion of using affine (linear) registration instead.

After calculating the transformation needed to register a given run to the reference run, we're ready to run the final transformation. This concatenates a number of transformations together, in order to keep the number of interpolations to a minimum. The steps are:

• motion correction
• unwarping
• registration to the reference run

Finally, the transformed functional timeseries is corrected for average bias field. This changes the mean intensity of the image at different voxels, but does not change the shape of the timeseries at any voxel.

Let's look at the commands that will be run to process one run:

reg_unwarp_bold_run.py bender_03 prex_1 study_1 -k --dry-run

As before, we'll set the -k flag to keep intermediate files, for illustrative purposes.

Note: depending on the details of your study, you may want to use just linear registration. To do this, add -l to the call to reg_unwarp_bold_run.py. That will cause registration to be affine (transformation + scaling).

To process all of the runs, for both subjects:

rlaunch "reg_unwarp_bold_run.py {s} {r} study_1 -k" bender_03:bender_04 prex_1:study_1 -N 1 -n 4 -a 6 -r 01:00:00 -p development

Note that we use the same script to process the study_1 run. Since it is the reference run, there will be no co-registration step in that case.

When it finishes, you should now have a new directory for each subject in BOLD/antsreg/data, which contains two images, prex_1.nii.gz and study_1.nii.gz. These are the functional timeseries with all basic preprocessing applied, in the space of the reference run. They are now ready for use with level 1 GLM modeling through FEAT.

Next Steps

See the FEAT documentation for details on univariate analysis in FSL. There are also scripts in FAT to help with running the various levels of univariate analysis, detailed on the main wiki.

Filtering and smoothing

If you're just running a GLM, you can easily apply filtering and smoothing when running your FEAT level 1 model. However, before running pattern classification analysis or representational similarity analysis, you may want to apply a temporal filter and smooth your BOLD data. See Filtering and Smoothing for details.