DTI Tutorial 2 - Normalization and Statistics
In our first tutorial we processed our raw images and calculated tensor images. In this second tutorial we will address how to normalize a set of diffusion tensor images (DTI) and then run statistics on the normalized brain images. First it is important to understand why you have to normalize brain images; lets see what Wikipedia has to say about spatial normalization:
"In neuroimaging, spatial normalization is an image processing step, more specifically an image registration method. Human brains differ in size and shape, and one goal of spatial normalization is to deform human brain scans so one location in one subject's brain scan corresponds to the same location in another subject's brain scan." - Wikipedia
As you see, inter-individual statistics will only be effective when brain images are in the same space, and spatial normalization can take care of that. An important distinction between diffusion tensor imaging and other imaging modality (like T1 anatomicals or fMRI scans), is that DTI allows for a more sophisticated normalization. Where T1 and fMRI scans are normalized based on intensity scales, diffusion tensor images instead offer the ability to register tensors, that have both length and direction. The property thus allows for a better capability to register white matter structures that have highly anisotropic tensors. In this video the developer of a publicly available tensor based normalization toolbox, Gary Zhang, discusses why tensor-based normalization can be better than other normalization methods, what statistical analyses you can run after normalization, and the difference between tract based spatial statistics (TBSS) and tract specific analysis (TSA). The below image details some of the advantages. For example in the area's with the yellow circles, you see areas of the brain where T1 anatomicals would have a hard time effectively registering inter-subjects anatomy for different white matter tracts, because the contrast in that region is equal despite it including multiple tracts with different orientation. Using the tensor information as show on the right will help discriminate between those regions, and allow for a more advanced registration between subjects.
In this tutorial there will be a number of steps that we have to run through in order to test statistical models, these steps are detailed below:
 |
| Advantage of tensor specific normalization. Image by DTI-TK |
1. Make sure the tensor files are in the right format for tensor normalization.
2. Make population template by normalizing tensor images
3. Register population template to MNI space
4. Produce scalar images (FA, MD, AD, RD) and check quality
5. Run whole brain voxel-wise statistics
6. Share your data.
4. Produce scalar images (FA, MD, AD, RD) and check quality
5. Run whole brain voxel-wise statistics
6. Share your data.
You can look through slides that go into some of the background of this tutorial here. Throughout the tutorial there will be questions to help with a better understanding of each of the steps and introduce you to the software that is being used. This might require installing said software and downloading example data from my website. At the end of the document you will find a pdf with the answers to the questions in this tutorial, although try to figure things out yourself as much as possible. If you run into issues going through this tutorial please direct your questions to the diffusion imaging forum, where there will likely be someone able to help you out. Please make sure to describe your question as clearly as possible. Remember to please direct software specific questions to their own respective forums or mail lists.
1. Standardize Tensor Images for Normalization
To get the tensor images ready for normalization we will have to make sure that tensors:
- have the right diffusivity unit,
- remove extreme outliers,
- and are symmetric and positive-definite.
- have the right diffusivity unit,
- remove extreme outliers,
- and are symmetric and positive-definite.
These steps are described in detail on the website of DTI-TK; which is the website of the software package that we will use for the tensor normalization.
Example code:
i=052212_s09_dti_eddy_fm_strip_dt.nii
prefix=052212_s09_dti_eddy_fm_strip
echo ~~~Adjusting Diffusivity Units~~~
TVtool -in $i -scale 1000000000 -out ${prefix}_sdt.nii.gz
echo ~~Checking for and removing Outliers~~
TVtool -in ${prefix}_sdt.nii.gz -norm
SVtool -in ${prefix}_sdt_norm.nii.gz -stats
BinaryThresholdImageFilter ${prefix}_sdt_norm.nii.gz ${prefix}_non_outliers.nii.gz 0 100 1 0
TVtool -in ${prefix}_sdt.nii.gz -mask ${prefix}_non_outliers.nii.gz -out ${prefix}_tmp.nii.gz
mv -f ${prefix}_tmp.nii.gz ${prefix}_sdt.nii.gz
TVtool -in ${prefix}_sdt.nii.gz -norm
echo ~Stats for${prefix} - max should be below 100~
SVtool -in ${prefix}_sdt_norm.nii.gz -stats
echo ~~~Enforcing positive semi-definiteness~~~
TVtool -in ${prefix}_sdt.nii.gz -spd -out ${prefix}_tmp.nii.gz
mv -f ${prefix}_tmp.nii.gz ${prefix}_sdt.nii.gz
echo ~~~Standardizing Voxel Space~~~
TVAdjustVoxelspace -in ${prefix}_sdt.nii.gz -origin 0 0 0 -out ${prefix}_sdt.nii.gz
echo ~~~Resample to isotropic voxels~~~
TVResample -in ${prefix}_sdt.nii.gz -align center -size 256 256 128 -vsize 1 1 1
rm -f ${prefix}_tmp.nii.gz
rm -f ${prefix}_non_outliers.nii.gz
rm -f ${prefix}_sdt_nonSPD.nii.gz
rm -f ${prefix}_sdt_norm.nii.gz
Q1: Go through the DTI-TK preprocessing tutorial, indicate which lines of the above code correspond to which steps in the outline. What is the goal of each step? When you use this code make sure you use the right diffusivity units, the code above is just an example, find out how to convert your own here.
This code uses a number of new tools. One of which will be of interest to take a closer look at.
TVtool by DTI-TK (download here)
http://dti-tk.sourceforge.net/pmwiki/pmwiki.php?n=Documentation.TVtool
Q2: Using TVtool investigate how and then write the code to produce the axial diffusivity (AD) map from the standardized tensor file (download 052212_s09_dti_eddy_fm_strip_sdt.nii.gz here). What other maps can that tool help produce?
Q3: What are the 4 main diffusion measurements that are often used in DTI? Take a look at the explanation on this page, and describe which measure you would pick if you were selectively interested in white matter?
This code uses a number of new tools. One of which will be of interest to take a closer look at.
TVtool by DTI-TK (download here)
http://dti-tk.sourceforge.net/pmwiki/pmwiki.php?n=Documentation.TVtool
Q2: Using TVtool investigate how and then write the code to produce the axial diffusivity (AD) map from the standardized tensor file (download 052212_s09_dti_eddy_fm_strip_sdt.nii.gz here). What other maps can that tool help produce?
Q3: What are the 4 main diffusion measurements that are often used in DTI? Take a look at the explanation on this page, and describe which measure you would pick if you were selectively interested in white matter?
Additional online resources:
- More info on how to calculate the diffusion tensor here
- Articles demonstrating improved registration when using DTI-TK here and here
- Articles demonstrating improved registration when using DTI-TK here and here
2. Normalize Tensors Images
In the next step we will take the standardized tensor images for each individual subject and run a multi-step normalization. The steps of the normalization are detailed on the DTI-TK website, but essentially consists of four steps:1. Make initial mean
2. Run a rigid registration
3. Run an affine registration
4. Run a diffeomorphic registration.
Each normalization step iteratively moves each individual image closer to the population mean. As described by DTI-TK:
"Rigid alignment aims to find a rigid transformation, a linear transformation that does not alter objects size and shape, that best aligns two objects in question. (...) Affine alignment improves upon rigid alignment by finding a more general linear transformation that allows objects' global size and shape to be altered. (...) After you've aligned your DTI volumes affinely, you can proceed with the deformable alignment, which is designed to improve alignment quality by removing size or shape differences between local structures." - DTI-TK
 |
Iterative normalization process. Image by DTI-TK
|
This process starts by averaging all the individual tensor images to create a starting point template, a population template that has yet to be normalized, the example FA from these tensor images may look something like this:
 |
| Example of un-normalized average FA of population of images |
Then as a next step each individual gets registered to this new population template iteratively. This process consists of rigid, affine and diffeomorphic registrations. A rigid-body transformation in three dimensions (x, y and z) is defined by six parameters: three translations and three rotations (source).
These transformations are a subset of the 12 transformations allowed in an affine registration of a three dimensional volume, which besides the rotation and translation also includes three sheers, and three zooms. Although affine registrations still have the constraint that parallel lines remain parallel after registration, thus it is contained to a linear registration. Or as wiki puts it:
 |
| Affine includes 12 parameter transformation: 3 rotations, 3 translations, 3 sheers, 3 zooms |
"In geometry, an affine transformation, affine map or an affinity (from the Latin, affinis, "connected with") is a function between affine spaces which preserves points, straight lines and planes. Also, sets of parallel lines remain parallel after an affine transformation. An affine transformation does not necessarily preserve angles between lines or distances between points, though it does preserve ratios of distances between points lying on a straight line." - Wikipedia
Finally the diffeomorphic or deformable registration removes those constraints too and allows for a non-linear higher dimension registration, as long as it is invertible.
Example script:
 |
| Examples of diffeomorphic transformations |
echo ~~~Create an initial mean~~~;
for spd in `ls *sdt.nii.gz`;
do
subj_list_file=subj_list;
echo $spd >> subj_list_file.txt;
done
TVMean -in $subj_list_file.txt -out mean_initial.nii.gz;
TVResample -in mean_initial.nii.gz -vsize 1 1 1 -size 256 256 128;
echo ~~~Rigid alignment with initial mean~~~;
sh dti_rigid_population mean_initial.nii.gz ${subj_list_file}.txt EDS 3;
echo ~~~Affine alignment with output from rigid alignment~~~;
sh dti_affine_population mean_rigid3.nii.gz ${subj_list_file}.txt EDS 3;
echo ~~~Diffeomorphic alignment with output from affine alignment~~~;
TVtool -tr -in mean_affine3.nii.gz;
BinaryThresholdImageFilter mean_affine3_tr.nii.gz mask.nii.gz 0 .01 100 1 0;
sh dti_diffeomorphic_population mean_affine3.nii.gz ${subj_list_file}_aff.txt mask.nii.gz 0.002
Q4: Which step do you think will take the most time to run? Note the step where you resample your image to have isotropic voxel size, which line is this? Read this post on the importance, explain why we need isotropic voxels for fiber tracking. Resampling also alters the image dimensions, when you change those, you can read here that DTI-TK needs your image dimensions to be powers of 2. What are most likely often used powers of 2 in neuroimaging volume dimensions?
Additional online resources:
- Courses on fMRI spatial normalization by Wellcome Trust centre here and here- Courses on fMRI spatial normalization by SPM here
- Advanced normalization tools (ANTS) for normalization of scalar images
- Background documentation on normalization by ANTS
- Further documentation from DTI-TK on tensor normalization
3. Register population template to MNI-atlas
Below you can see the result of the population mean post-normalization. At this stage you are at a cross-roads of options how to analyze the results. The following steps are but one way, as you get more comfortable with analyses, software and scripting you can try to explore other options.
 |
| Example of normalized average FA of population of images |
 |
| MNI ICBM152 non-linear 6th generation. Image by FSL. |
 |
| Template FA (red) after registration to MNI space (gray-scale) |
Example script:
echo ~~~Register population template to MNI~~~
TVtool -in ${population_tensor_mean} -tr -out my_tr.nii.gz
BinaryThresholdImageFilter my_tr.nii.gz my_mask.nii.gz 0.7 100 1 0
TVtool -in ${population_tensor_mean} -mask my_mask.nii.gz -out ${population_tensor_mean}_strip.nii.gz
TVtool -in ${population_tensor_mean}_strip.nii.gz -fa -out ${population_tensor_mean}_strip_fa.nii.gz
asvDSM -template ${MNI_T1} -subject ${population_tensor_mean}_strip_fa.nii.gz -outTrans ${dtiToT1Trans} -sep 0.5 0.5 0.5 -ftol 0.0001
echo ~~~Run subject specific loop to warp to T1~~~
for i in `cat subjects.txt`;
do
echo ~Compose full warp and apply~
dfRightComposeAffine -aff ${subj}_affineTrans -df ${subj}_diffeoTrans -out ${subj}_combined.df.nii.gz
echo ~Combine full warp with affine to T1 and apply~
dfLeftComposeAffine -df ${subj}_combined.df.nii.gz -aff ${dtiToT1Trans} -out ${subj}_toT1_combined.df.nii.gz
deformationSymTensor3DVolume -in ${subj}_orig -trans ${subj}_toT1_combined.df.nii.gz -target ${T1} -out ${subj}_tensorFile
done
Q5: Open fslview, click "File" -> "Open standard". In this directory and in the directory atlases you can browse all the different atlases that are automatically loaded with FSL. You can also find the atlas images in your installed FSL directory. Find the MNI template that you would like to use for registration. What is the name of that file? Give 1 reason why you would want to register your data to this template.
Additional online resources:
- More on the available atlases in FSL- More on the atlases by MNI
- White matter atlas by Hermoye et al.
- International consortium of brain mapping (ICBM) DTI-81 atlas
- Cortical association tracts by Mori et al. (2002)
- Fiber tract based atlas of human white matter by Wakana et al. (2004)
- A diffusion tensor tractography atlas by Catani & De Schotten (2008)
- More online atlases by JHU
4. Make scalars and check quality
As a next step we will make diffusion maps (FA, MD, AD, RD) of each individual in template space. These diffusion maps are also known as scalar images; they have magnitude but not direction. When you have those images you can go through them and get an impression of the quality of the normalization and your data. One way of quickly viewing all your images is by concatenating them into one 4 dimensional image using fslmerge.
fslmerge by FSL
http://fsl.fmrib.ox.ac.uk/fsl/fslwiki/Fslutils
Example script:fslmerge by FSL
http://fsl.fmrib.ox.ac.uk/fsl/fslwiki/Fslutils
 |
| Movie of normalized FA images |
echo ~~~Produce subject specific smoothed scalars~~~
for subj in `cat subjects.txt`;
do
echo ~Produce diffusion maps~
TVtool -in ${subj}_tensorFile.nii.gz -fa -out ${subj}_fa.nii.gz
TVtool -in ${subj}_tensorFile.nii.gz -tr -out ${subj}_tr.nii.gz
TVtool -in ${subj}_tensorFile.nii.gz -ad -out ${subj}_ad.nii.gz
TVtool -in ${subj}_tensorFile.nii.gz -rd -out ${subj}_rd.nii.gz
echo ~Smooth diffusion maps~
SVGaussianSmoothing -in ${subj}_fa.nii.gz -fwhm 4 4 4 -out ${subj}_fa4mm.nii.gz
SVGaussianSmoothing -in ${subj}_tr.nii.gz -fwhm 4 4 4 -out ${subj}_tr4mm.nii.gz
SVGaussianSmoothing -in ${subj}_ad.nii.gz -fwhm 4 4 4 -out ${subj}_ad4mm.nii.gz
SVGaussianSmoothing -in ${subj}_rd.nii.gz -fwhm 4 4 4 -out ${subj}_rd4mm.nii.gz
done
echo ~~~Merge FA images for quality control~~~
fslmerge -t all_images_MNI_FA *_fa.nii.gz
Q6: Some of the code in this tutorial uses a "for loop" for the first time. Take a look at a tutorial on loops here. Try to figure out how "for loops" work, then describe what you would put in the subjects.txt file from line 2. HINT: Remember how variables work in bash. Another option is to run the below code, instead of line 2 and 3. Can you describe what this would do?
Q8: Next describe what lines 10-13 accomplish, if you are unsure read this documentation on spatial smoothing. Finally in line 16 we merge all the individual FA images, describe why you might want to do that.
for i in `ls *_diffeo_fa.nii.gz`;
do
subj=`echo $i | awk 'BEGIN{FS="_diffeo_fa.nii.gz"}{print $1}'`
Q7: Note that DTI-TK outputs the trace (TR) instead of the MD value. Look up how trace relates to MD here, and write the fslmaths code to convert TR to MD. If you are unsure how to do this, type fslmaths in your command line to see all the help options. Q8: Next describe what lines 10-13 accomplish, if you are unsure read this documentation on spatial smoothing. Finally in line 16 we merge all the individual FA images, describe why you might want to do that.
Additional online resources:
- More info on smoothing your data- How to convert FWHM to sigma
5. Run whole brain voxel-wise statistics
Now that we have scalar images in normalized space for each individual we can run voxel-wise analyses. There is a plethora of tools that you can use to run your preferred statistical model. A method that is often used in diffusion imaging is tract-based spatial statistics (TBSS) by FSL. Although this method has become the somewhat accidental standard, it has recently come under some scrutiny. Bach et al. (2014) describe some important methodological considerations for using TBSS. And Swartz et al. (2014) describe how TBSS's skeleton projection step actually reduces algorithm accuracy. This is due to developments in registration tools (e.g. tensor normalization with DTI-TK, and improved algorithms of advanced normalization tools [ANTS]), which has made the methods used in TBSS obsolete. Furthermore, there is the possibility that you miss subtle effects that are in one of the smaller tracts, that may get removed when the data is skeletonized.
Due to these issues, I have a preference for running whole brain voxel-wise analyses. In this case I decided to use FSL's randomise, since it is free and does not depend on software that has to run in MatLab (like SPM and SurfStat). Randomise runs a non-parametric test, which means that it makes no assumptions about the underlying normality of the data, which is a more robust test than ordinary GLM's. If you read the randomise documentation you will also see that you can include covariates like age and sex. Which are important variables to include in anatomical studies.
 |
| Example of tract-based spatial statistics. Image by FSL |
randomise by FSL
http://fsl.fmrib.ox.ac.uk/fsl/fslwiki/Randomise/UserGuide
Example code:
randomise -i all_merge_fa.nii.gz -o randomise_out -d test.mat -t test.con -m my_mask.nii.gz -T
Q9: Look at the user guide for randomise, and try to get an impression of what the .mat file and .con file contain. Produce an example of both for the imaginary case that you have a patient group (n=3) and healthy control group (n=3), and you want to test where FA is significantly lower in the patient group, and where FA is significantly lower in the control group. What does the option -T stand for in randomise?
Q10: Open fslview and open a standard 1mm MNI T1 image. Download example t-stat map (test_tstat1.nii.gz here), and add it to the MNI image. What you see is not yet very informative, click on the information dialog button (i in blue circle) and change the Lookup table options. Make the first Red-Yellow, and turn the second one on, and make it Blue-Lightblue. Close the dialog screen. Go to this p-value calculator, and input in the t-ratio boxes that you want a probability that is 0.05 and you have 5 degrees of freedom (number of subjects - 1). What is your 2-tailed significant t-statistic? Input that number in the fslview screen as a minimum brightness value, and make 5 your maximum brightness value. Although this is obviously just test data, you now see what a noisy (when you only test 6 subjects you will most likely only find noise) t-stat map can look like. It is now up to you to start using your own data and improve on this! Remember to correct for multiple comparison when running your data to avoid inflated false discovery rates. You are running a t-tests on each individual voxel, which combine to multiple thousands per brain, and you will have to correct for this. Read more about familywise error rates here. What does the -T option do to help with this issue? HINT: Read this .pdf on the background of threshold-free cluster enhancement (TFCE).
Q10: Open fslview and open a standard 1mm MNI T1 image. Download example t-stat map (test_tstat1.nii.gz here), and add it to the MNI image. What you see is not yet very informative, click on the information dialog button (i in blue circle) and change the Lookup table options. Make the first Red-Yellow, and turn the second one on, and make it Blue-Lightblue. Close the dialog screen. Go to this p-value calculator, and input in the t-ratio boxes that you want a probability that is 0.05 and you have 5 degrees of freedom (number of subjects - 1). What is your 2-tailed significant t-statistic? Input that number in the fslview screen as a minimum brightness value, and make 5 your maximum brightness value. Although this is obviously just test data, you now see what a noisy (when you only test 6 subjects you will most likely only find noise) t-stat map can look like. It is now up to you to start using your own data and improve on this! Remember to correct for multiple comparison when running your data to avoid inflated false discovery rates. You are running a t-tests on each individual voxel, which combine to multiple thousands per brain, and you will have to correct for this. Read more about familywise error rates here. What does the -T option do to help with this issue? HINT: Read this .pdf on the background of threshold-free cluster enhancement (TFCE).
Additional online resources:
- YouTube channel on analysis of neuroimaging data by Jeanette Mumford- Tutorial on including voxel-wise covariates using SurfStat
- Stanford university information on voxel-wise tensor statistics and whole brain diffusion stats
- Wikipedia on false discovery rate (FDR)
- Lecture on Multiple comparison issues from Columbia University
- Genovese et al. (2002) on FDR for neuroimaging
- Brainvoyager software on multiple comparison correction and FDR
- More info on thresholding by the late Keith Worthsley, and the seminal publication on significance testing in MRI by Worthsley et al. (1996)
- Wikipedia on false discovery rate (FDR)
- Lecture on Multiple comparison issues from Columbia University
- Genovese et al. (2002) on FDR for neuroimaging
- Brainvoyager software on multiple comparison correction and FDR
- More info on thresholding by the late Keith Worthsley, and the seminal publication on significance testing in MRI by Worthsley et al. (1996)
6. Share your data
When you are done with your analyses and want to share your statistical maps publicly, you can use NeuroVault. Sharing your statistical maps will help improve meta analyses, collaborations and science in general.
 |
| Upload your unthresholded statistical map here, and share your URL in your publication |
Answers:
DTI Tutorial 2 Answers---------------------
Reference this post as: Do Tromp, DTI Tutorial 2 - Normalization and Statistics, The Winnower 3:e145653.31599 (2016). DOI: 10.15200/winn.145653.31599
