If You're at an Organized Course
If you are taking one of the formally organized courses, everything has been set up for you on the provided laptop. The only thing you will need to do is run the following commands in every new terminal window (aka shell) you open throughout this tutorial. Copy and paste the commands below to get started:
export SUBJECTS_DIR=$TUTORIAL_DATA/buckner_data/tutorial_subjs/group_analysis_tutorial cd $SUBJECTS_DIR/glm
To copy: Highlight the command in the box above, right click and select copy (or use keyboard shortcut Ctrl+c), then use the middle button of your mouse to click inside the terminal window to paste the command (or use keyboard shortcut Ctrl+Shift+v). Press enter to run the command.
These two commands set the SUBJECTS_DIR variable to the directory where the data is stored and then navigates into the sub-directory you'll be working in. You can now skip ahead to the tutorial (below the gray line).
If You're Not at an Organized Course
If you are NOT taking one of the formally organized courses, then to follow this exercise exactly be sure you've downloaded the tutorial data set before you begin. If you choose not to download the data set, you can follow these instructions on your own data, but you will have to substitute your own specific paths and subject names. These are the commands that you need to run before getting started:
## bash <source_freesurfer> export TUTORIAL_DATA=<path_to_your_tutorial_data> export SUBJECTS_DIR=$TUTORIAL_DATA/buckner_data/tutorial_subjs/group_analysis_tutorial cd $SUBJECTS_DIR/glm ## tcsh source $FREESURFER_HOME/SetUpFreeSurfer.csh setenv TUTORIAL_DATA <path_to_your_tutorial_data> setenv SUBJECTS_DIR $TUTORIAL_DATA/buckner_data/tutorial_subjs/group_analysis_tutorial cd $SUBJECTS_DIR/glm
If you are not using the tutorial data, you should set your SUBJECTS_DIR to the directory in which the recon(s) of the subject(s) you will use for this tutorial are located.
This tutorial is designed to introduce you to the "command-line" group analysis stream in FreeSurfer (as opposed to QDEC which is GUI-driven). While this tutorial shows you how to perform a surface-based thickness study, it is important to realize that most of the concepts learned here apply to any group analysis in FreeSurfer, surface, volume, thickness or fMRI. Here are some useful Group Analysis links you might want to refer back to at a later time:
DODS vs DOSS
This Data Set
The data used for this tutorial is 40 subjects from Randy Buckner's lab. It consists of males and females, ages 18 to 93. You can see the demographics here. You will perform an analysis looking for the effect of age on cortical thickness, accounting for the effects of gender in the analysis. This is the same data set used in the QDEC Tutorial, and you will get the same result.
General Linear Model (GLM) DODS Setup
Design Matrix/FSGD File
In this tutorial, we will model the thickness as a straight line. A line has two parameters: an intercept (or offset) and a slope. For this example:
- The slope is the change of thickness with age.
- The intercept/offset is interpreted as the thickness at age=0.
Parameter estimates are also called "regression coefficients" or "betas". To account for effects of gender, we will model each sex with its own line, meaning that there will be four linear parameters:
- Intercept for Females
- Intercept for Males
- Slope for Females
- Slope for Males
In FreeSurfer, this type of design is called DODS (for "Different-Offset, Different-Slope").
You can either create your own design matrices, or, if you specify your design as a FreeSurfer Group Descriptor File (FSGD), FreeSurfer will create the design matrices for you. The FSGD file is a simple text file you create. It is not generated by FreeSurfer. See this page for the format. The demographics page also has an example FSGD file for this data.
Exercises to try on your own:
- Create an FSGD file for the above design. For the tutorial, one (gender_age.fsgd) already exists so that you can continue with the exercises.
Information on how to create/view text files can be found here.
A contrast is a vector that embodies the hypothesis we want to test. In this case, we wish to test the change in thickness with age, after removing the effects of gender. To do this, create a simple text file with the following numbers (if you're using the tutorial data, this has already been done for you):
0 0 0.5 0.5
- Remember that for an analysis between two groups each with an intercept of b and a slope of m, the contrast matrix will be in the format [b1 b2 m1 m2].
- There is one value for each parameter (so 4 values total).
- The intercept/offset values (b1, b2) are 0 (nuisance).
- The slope values (m1, m2) are 0.5 so as to average the Female and Male slopes.
You'll find we created the contrast matrix for you already for this tutorial. It's a file called lh-Avg-thickness-age-Cor.mtx. You will need to create your own contrast matrix when testing your own hypotheses on your data. FreeSurfer doesn't automatically create this file. There are examples of FSGD files and contrast matrices for many different hypotheses here.
Assemble the Data (mris_preproc)
Assembling the data simply means:
- Resampling each subject's data into a common space.
- Concatenating all the subjects' into a single file.
- Spatial smoothing (can be done between 1 and 2).
This can be done in two equivalent ways:
1. Previously Cached (qcached) Data
When you run recon-all with the -qcache option (see here for more info), recon-all will resample data onto the average subject (called fsaverage) and smooth it at various FWHM (full-width/half-max) values, usually 0, 5, 10, 15, 20, and 25mm. This can speed up later processing. The data for this tutorial was run with the -qcache option, so the FWHM are already cached. Run this command to assemble your previously cached data:
mris_preproc --fsgd gender_age.fsgd \ --cache-in thickness.fwhm10.fsaverage \ --target fsaverage \ --hemi lh \ --out lh.gender_age.thickness.10.mgh
NOTE: The backslash allows you to copy and paste multiple lines of code as one command. We use this throughout the tutorials to display the commands in a more easy-to-read manner, while still allowing you to copy and paste. Whenever you are typing in your own commands, instead of copying and pasting a command written out across multiple lines, a backslash is not necessary.
- Only takes a few seconds because the data have been cached.
- The FSGD file lists all the subjects, helps keep order.
- The independent variable is the thickness smoothed to 10mm FWHM.
The target is the average subject, fsaverage. This subject is included with the FreeSurfer distribution.
- Analyses will only be run on the left hemisphere.
Since the data is already cached, this mris_preproc command just creates a stack of all the smoothed thickness data from each subject into one file as output.
The output of this command is given by the --out flag. Here we call it lh.gender_age.thickness.10.mgh.
Things to do:
- Run mri_info on lh.gender_age.thickness.10.mgh to find its dimensions:
2. Uncached Data
In the case that you have not cached the data, you can perform the same operations manually using the two commands below. These two commands give the equivalent results as running recon-all with -qcache. However, the commands below only create output for unsmoothed data and data smoothed to 10mm FWHM. qcache, on the other hand, smooths the data for a range of different FWHM vaules.
THESE NEXT TWO COMMANDS ARE OPTIONAL FOR THE TUTORIAL, AND WILL TAKE ABOUT 5 MINUTES TO RUN.
The tutorial data has already been pre-processed (cached) for you.
mris_preproc --fsgd gender_age.fsgd \ --target fsaverage \ --hemi lh \ --meas thickness \ --out lh.gender_age.thickness.00.mgh
- This resamples each subject's left hemisphere data to fsaverage.
- Output is lh.gender_age.thickness.00.mgh, which is unsmoothed.
mri_surf2surf --hemi lh \ --s fsaverage \ --sval lh.gender_age.thickness.00.mgh \ --fwhm 10 \ --cortex \ --tval lh.gender_age.thickness.10B.mgh
- This smooths each subject's resampled data by 10mm FWHM.
- "--cortex" means only smooth areas in cortex (exclude medial wall). This is automatically done with qcache. You can also specify other labels.
- "--sval" is the name of the source subject found in the $SUBJECTS_DIR. This input data must be sampled onto this subject's surface.
- "--tval" is the name of the file where the data on the target surface will be stored.
- Output is lh.gender_age.thickness.10B.mgh. (We have named the output file ...10B.mgh because the precreated file from running qcache ...10.mgh already exists.)
GLM Analysis (mri_glmfit)
mri_glmfit \ --y lh.gender_age.thickness.10.mgh \ --fsgd gender_age.fsgd dods\ --C lh-Avg-thickness-age-Cor.mtx \ --surf fsaverage lh \ --cortex \ --glmdir lh.gender_age.glmdir
- Input is lh.gender_age.thickness.10.mgh.
- Same FSGD used as with mris_preproc. Maintains subject order!
- DODS is specified (it is the default).
- Only one contrast is used (lh-Avg-thickness-age-Cor.mtx), but you can specify multiple contrasts.
- "--cortex" specifies that the analysis only be done in cortex (ie, medial wall is zeroed out). Other labels can be used.
- The output directory is lh.gender_age.glmdir.
- Should only take about 1min to run.
Things to do:
When this command is finished the glm directory will contain a sub-directory called lh.gender_age.glmdir. There will be a number of output files in this directory, as well as other subdirectories. Run:
and you will see the following files (descriptions are also included below, but don't show up in your terminal screen):
beta.mgh -- all parameter estimates (surface overlay) dof.dat -- degrees of freedom (text) fwhm.dat -- average FWHM of residual (text) lh-Avg-thickness-age-Cor -- contrast subdirectory mask.mgh -- binary mask (surface overlay) mri_glmfit.log -- log file (text, send this with bug reports) rstd.mgh -- residual standard deviation (surface overlay) rvar.mgh -- residual variance (surface overlay) sar1.mgh -- residual spatial AR1 (surface overlay) surface -- the subject and hemisphere used for this analysis (text) Xg.dat -- design matrix (text) X.mat -- design matrix (MATLAB format) y.fsgd -- copy of input FSGD file (text)
NOTE: You may also see some temporary directories (with names starting with tmp) or log files, which you can ignore.
- The DOF is the degrees of freedom for the analysis (has important implications for the p-value).
- The FWHM is a measure of the smoothness of the data. Part of this comes from the applied smoothing (eg, 5mm FWHM) and part of the smoothness comes from the inherent smoothness in the data. When correcting for multiple comparisons, the final FWHM must be taken into account (but this is done automatically).
There will be a subdirectory for each contrast that you specify. The name of the directory will be that of the contrast matrix file (without the .mtx extension). For example, to inspect the lh-Avg-thickness-age-Cor directory, run
and you will see the following files:
C.dat -- original contrast matrix (text) cnr.mgh -- contrast-to-noise ratio (surface overlay) efficiency.dat -- statistical efficiency for the contrast (text) F.mgh -- F ratio of contrast (surface overlay) gamma.mgh -- contrast effect size (surface overlay) gammavar.mgh --contrast variance (surface overlay) maxvox.dat -- voxel with the maximum statistic (text) pcc.mgh -- partial (pearson) correlation coefficient (surface overlay) sig.mgh -- significance, -log10(pvalue), uncorrected (surface overlay) z.mgh -- z-stat that corresponds to the significance (surface overlay)
View the uncorrected significance map with freeview. First, make sure you are in the correct directory:
Then, run this command to visualize the data:
freeview -f $SUBJECTS_DIR/fsaverage/surf/lh.inflated:annot=aparc.annot:annot_outline=1:overlay=lh.gender_age.glmdir/lh-Avg-thickness-age-Cor/sig.mgh:overlay_threshold=4,5 -viewport 3d -layout 1
You should see something similar to this:
The threshold is set to 4, meaning vertices with p<.0001, uncorrected, will have color.
- The lower threshold sets the minimum sigificance that a voxel must meet in order to be visible. It should be set such that only likely true effects will be seen. You can kind of think of it as a filter that removes voxels where no effect is present. The upper is just used for visualization purposes. If you set it very close to the min, then all voxels will appear to be the same color. If you set it higher, then you can get some gradation due to the strength of the effect. So the min removes voxels without an effect and the max allows you to see how the strength of the effect varies across the brain.
- Blue represents a negative correlation (i.e. thickness decreases with age), red represents positive correlation (i.e. thickness increases with age).
- In this example, the sig.mgh is used as the overlay and in practice this is what is normally used. However, the pcc.mgh can also be used as the overlay if it was created in the analysis.
- Click on a point in the Precentral Gyrus. What is its value? What does it mean?
Viewing the medial surface, change the overlay threshold to something very, very low (say, .01), by clicking Configure (under the Overlay dropdown menu that says sig.mgh) and set the Min value to 0.01. Click Apply, or check the Automatically apply changes checkbox.
You should see something like this:
- Almost all of the cortex now has color.
- The non-cortical areas are still blank (0) because they were excluded with --cortex in mri_glmfit above.
All the surface overlays created by mri_glmfit, not just the significance map, can also be inspected in freeview. Simply replace the path to the sig.mgh file in the command above to the surface overlay that you would like to see. For example, this will display the F ratio instead of the significance map:
freeview -f $SUBJECTS_DIR/fsaverage/surf/lh.inflated:annot=aparc.annot:annot_outline=1:overlay=lh.gender_age.glmdir/lh-Avg-thickness-age-Cor/F.mgh:overlay_threshold=20,50 -viewport 3d -layout 1
Text output and log files can be opened in any text editor.
What is a contrast? Answer
Which type of file does FreeSurfer not automatically create? Why not?Answer
What colors are associated with a negative and positive correlation? Answer
What does the threshold of 4 mean? What would the threshold be if the p-value was < 0.001? What about < 0.05? Answer
What was the DOF for this experiment? Answer
What was the FWHM? Answer
How many frames does beta have? Why? (hint: run mri_info on lh.gender_age.glmdir/beta.mgh) Answer