Differences between revisions 1 and 5 (spanning 4 versions)
Revision 1 as of 2011-03-08 14:10:46
Size: 48741
Editor: tanha
Comment:
Revision 5 as of 2011-03-08 19:35:54
Size: 15582
Editor: tanha
Comment:
Deletions are marked like this. Additions are marked like this.
Line 22: Line 22:
 * Each subject has 1 run (except sess01 which has 4 runs)
Line 23: Line 24:
 * TR = 3 sec, 85 time points/run (run length = 255 sec)
Line 32: Line 32:
environment variable. NOTE: if you are taking the Multi-Modal
Imaging C
lass at MGH, the data have already been set up.

Cd into the tutorial data directory and set SUBJECTS_DIR to use this data:
environment variable. NOTE: if you are taking a class at MGH, the data
have already been set up on your computer.

cd into the tutorial data directory and run ls:
Line 39: Line 39:
setenv SUBJECTS_DIR $PWD
}}}

It will be assumed that you are in this directory (or subdirectory
thereof) throughout the tutorial. This directory will have five
folders:

 * fb1-raw - raw fMRI data
 * fb1-raw-study - raw data organized as an FSFAST study but unanalyzed
 * fb1-preproc-study - The same data preprocessed.
 * fb1-analysis-study - The same data analyzed
 * subjects - FreeSurfer reconstruction of anatomical data (plus the fsaverage subject).

If you look in fb1-raw you will see 20 NIFTI data sets:
{{{
ls fb1-raw
f.mgh-101.1-r003.nii f.mgh-103.1-r003.nii f.mgh-104.1-r003.nii f.mgh-105.1-r003.nii f.mgh-106.1-r003.nii
f.mgh-101.1-r005.nii f.mgh-103.1-r005.nii f.mgh-104.1-r005.nii f.mgh-105.1-r005.nii f.mgh-106.1-r005.nii
f.mgh-101.1-r007.nii f.mgh-103.1-r007.nii f.mgh-104.1-r007.nii f.mgh-105.1-r007.nii f.mgh-106.1-r007.nii
f.mgh-101.1-r010.nii f.mgh-103.1-r010.nii f.mgh-104.1-r010.nii f.mgh-105.1-r010.nii f.mgh-106.1-r010.nii
}}}

These are the 20 fMRI runs mentioned above. In the name,
f.mgh-10X.1-rYYY.nii, X indicates the subject number (1, 3, 4, 5, 6),
and YYY indicates the run number. The sensory-motor task happened to be
runs 3, 5, 7, and 10.

---------------------------------------------------------------------
---------------------------------------------------------------------
= Quick Visualization Tutorial (tkmedit/tksurfer) =

The purpose of this tutorial is to familiarize you with how to use
FreeSurfer volume viewer (tkmedit) and surface viewer (tksurfer) in
the context of viewing functional data. You should otherwise already
know how to use tkmedit and tksurfer. See the pages below for a more
detailed handling of tkmedit and tksurfer.

== More info on tkmedit ==
 * [[TkMeditGuide/TkMeditGeneralUsage/TkMeditInterface|Interface]] <<BR>>
 * [[FsTutorial/TkmeditGeneralUsage|General usage]] <<BR>>
 * [[FsTutorial/TkmeditWorkingWithData|Working with data]] <<BR>>
 * [[FsTutorial/TkmeditReference|Quick reference]] <<BR>>

== More info on tksurfer ==

 * [[TkSurferGuide/TkSurferGeneralUsage/TkSurferInterface|Interface]] <<BR>>
 * [[FsTutorial/TksurferGeneralUsage|General usage]] <<BR>>
 * [[FsTutorial/TksurferWorkingWithData|Working with data]] <<BR>>
 * [[FsTutorial/TksurferReference|Quick reference]] <<BR>>

== Viewing a single functional overlay in the volume ==

cd into the study with already analyzed data:
{{{
cd $FSFTUTDIR/fb1-analysis-study
}}}
Run tkmedit-sess (this is an FSFAST wrapper for tkmedit):
{{{
tkmedit-sess -s mgh-101.1 -a sm-gamma-fwhm5 -c odd-v-0 -aparc+aseg
}}}

Don't worry about what all the arguments mean, as this part is only about
visualization. This command will bring up two windows: one with a brain
image, the other a control panel.

{{attachment:tkm-101-gam-cor.2-4.gif}}
{{attachment:tkm-101-toolswindow.gif}}

 * The brain image is the FreeSurfer anatomical for this subject.
 * The slightly pale colors on the anatomical indicate the FreeSurfer automatic segmentation.
 * The bright red/yellow/blue are super-threshold voxels in the functional overlay.
 * As you click on different points, you will see the "Functional value" field in the Tools window change as well as the "Sgmtn label".
 * The interpretation of the value depends on what is being viewed. This is a significance map, so the value is -log10(pvalue)*sign (i.e., for a pvalue = .01, the functional value will be +2).
 * The sign is a functional direction. The red/yellow are positive, and the blue are negative. As a functional value gets more positive, its color will change from red to yellow. As it gets more negative, it will change from blue to cyan.
 * Note that areas that are not above threshold will still have non-zero functional values.

{{{
Toggle the functional overlay on and off by hitting the button with
the red/yellow blob in the Tools window (it's in the top row - if you
mouse over it, you'll see a "Show Functional Overlay" tooltip).
}}}

{{{
Configure the functional overlay by clicking on "View..." in the Tools
window, then "Configure->Functional Overlay...". You should see the
following interface:
}}}

{{attachment:tkm-101-config-funcoverlay.gif}}

 * The thresholds are currently set at 2 (Min) and 4 (max).
 * The Min threshold is the minimum absolute value needed for a voxel to show an overlay color.
 * The maximum is the value beyond which the voxel will stop changing color.

{{{
Try changing these values, then hit the Apply button to see their effect.
}}}

So, why does it look like there's activation OUTSIDE of the brain? Try hitting View->Mask Functional Overlay to Aux button to see a marked improvement. And why is it so blocky? There are big voxels in
functional data, but if you'd like to pretend otherwise, try hitting
the "Trilinear" button on the Configure Functional Overlay window.

== Viewing multiple functional overlays and time courses in the volume ==

Exit from your current tkmedit window. Run tkmedit-sess (this is an FSFAST wrapper for tkmedit):

{{{
cd $FSFTUTDIR
setenv SUBJECTS_DIR $PWD
cd $FSFTUTDIR/fb1-analysis-study
tkmedit-sess -s mgh-101.1 -a sm-fir-fwhm5 -c odd-v-0 -aparc+aseg
}}}

This will bring up the two windows as you saw before (the brain image
window will have a different overlay). It will also bring up a window
called "Time Course". Click anywhere in the volume, and you will see a
time course associated. If you happened to have clicked in some activation, then
the time course will look similar to the one below.

{{attachment:tkm-101-fir.gif}}

The meaning of the time course depends upon the nature of the time
course loaded. In this case it is the hemodynamic response to the task
block averaged over all blocks over all runs (this is an "FIR"
model). The horizontal axis is time (0 means the onset of the
block). You can visualize the values of any multi-frame data set in
this way (it does not have to be a "time" course).<<BR>>

Notice that there is not much activation in the functional
overlay. This is good! You are looking at an overlay that corresponds
3 seconds prior to stimulus onset, so there should be no
activation.
{{{
To see the other time points, bring up the functional overlay configuration window
(View->Configure->Functional Overlay). Notice that the "Time Point" field now has a
range of 0-8 indicating the 9 time points in the Time Course window. If you hit the
+ button next to "Time Point" then hit Apply, you will see the overlay change.
}}}

You will also see a vertical dashed line in the Time Course
window. You are now looking at the map associated with the time between 0
and 3 sec after stimulus onset.

{{{
Keep hitting the +/Apply buttons to see different time points. You can
hit the "|>" button to start a movie of the activation (then "[]" to stop it).
}}}

== Viewing a single functional overlay on the surface ==

To view the same data on the surface, start a new sheel and run:

{{{
cd $FSFTUTDIR/fb1-analysis-study
tksurfer-sess -s mgh-101.1 -a sm-gamma-fwhm5 -c odd-v-0 -hemi lh -aparc
}}}

Again you will see two windows, "Tksurfer Tools" and a surface image
window like the one below. As with tkmedit, you can configure the
overlay with View->Configure...->Overlay. The second command will load
the surface parcellation which will hide the overlay. Change it to
outline mode with: View->Label Style->Outline.

{{attachment:tks-101-gam-lat.2-4.gif}}

You can view the time course with:
{{{
tksurfer-sess -s mgh-101.1 -a sm-fir-fwhm5 -c odd-v-0 -hemi lh -aparc
}}}

-------------------------------------------------------------------
-------------------------------------------------------------------
= Assembling the Data into the FSFAST Hierarchy =

Before analyzing data in FSFAST, it needs to be in a "Hierarchy", i.e., a
certain directory structure and naming convention. By adhering to a
hierarchy, you allow FSFAST to analyze large amounts of data
automatically because it knows where to find the input and where to
put the output.

 * Each step in the analysis modifies the hierarchy in some way.
 * All the data collected in one functional session (or visit) are placed under the folder (the "session directory").
 * All the sessions are arranged under one folder (called the "Study Directory").

The full hierarchy is shown in the image below:

{{attachment:fsfast-hierarchy.jpg}}

The red boxes indicate a directory with raw or preprocssed functional
data. E.g., f.nii is a raw fMRI 4D file, and fmc.nii is motion
corrected.

== Creating a Hierarchy from the Tutorial Data ==

There are five subjects of data, each with four runs. We will create a
an FSFAST session from one of the subjects. A properly created
hierarchy can be found in $FSFTUTDIR/fb1-raw-study.

Create a Study Directory and cd into it:
{{{
cd $FSFTUTDIR
mkdir -p my-study
cd my-study
}}}

Create a session directory for the first subject:
{{{
mkdir -p mgh-101.1
}}}

Create a bold and four run directories:
{{{
mkdir -p mgh-101.1/bold
mkdir -p mgh-101.1/bold/003
mkdir -p mgh-101.1/bold/005
mkdir -p mgh-101.1/bold/007
mkdir -p mgh-101.1/bold/010
}}}

Copy the raw data:
{{{
cp ../fb1-raw/f.mgh-101.1-r003.nii mgh-101.1/bold/003/f.nii
cp ../fb1-raw/f.mgh-101.1-r005.nii mgh-101.1/bold/005/f.nii
cp ../fb1-raw/f.mgh-101.1-r007.nii mgh-101.1/bold/007/f.nii
cp ../fb1-raw/f.mgh-101.1-r010.nii mgh-101.1/bold/010/f.nii
}}}
Note that they are all named "f.nii" but are distinguished by the
directory they are in.

== Link to FreeSurfer Anatomical Analysis ==

FSFAST makes extensive use of the FreeSurfer anatomical analysis. To
link a subject's functional and structural data, create a "subjectname"
file in the session directory. This file is literally called
"subjectname" and its contents have the FreeSurfer subject name
assigned during reconstruction and found in $SUBJECTS_DIR. These
subjects have already been reconstructed; this subject's name is
"fbph1-101". You can create this file with any text editor (emacs, vi,
pico), or just:

{{{
cd $FSFTUTDIR/my-study
echo fbph1-101 > mgh-101.1/subjectname
}}}

== Create a Stimulus Schedule (Paradigm File) ==

A "paradigm" file is a record of which stimulus was presented when and
for how long. Each paradigm file has four columns (extras are
ignored). The columns are:

 * 1. Stimulus onset time (sec)
 * 2. Condition ID code (0, 1, 2, ...)
 * 3. Stimulus Duration (sec)
 * 4. Stimulus Weight (usually 1)

The Stimulus onset is with respect to the time of the first image. The
Condition is a unique number that identifies the event or block type
that the stimulus belongs to. They must be contiguous, and Condition 0
is reserved for the Null Event (usually Fixation). The weight can be
used to model various effects of learning or adaptation (just set them
to 1 for now).

For this exercise, you will create two paradigm files based on the
fBIRN sensory-motor experiment. The timing for the experiment was
described above. The task timing is shown graphically below.
 
{{attachment:fbirn1-taskplot.gif}}

{{{
Using a text editor (emacs, vi, pico, etc), create a paradigm file called "sensory-motor0.par". Eg:
cd $FSFTUTDIR/fb1-raw-study
emacs sensory-motor0.par&
}}}
The correct paradigm file is here SensMotPar (This example file has a
5th column with the condition name. This is ignored in the analysis).

{{{
To make things a little more interesting, we are going to create
another paradigm file pretending that the odd and even blocks are
different. Assign Odd blocks Condition 1 and Even blocks Condition 2.
cd $FSFTUTDIR/fb1-raw-study
emacs sensory-motor.par&
}}}
Call this paradigm file "sensory-motor.par". The correct paradigm file
here is SensMotOddEvenPar.

Now copy this file into each run directory:
{{{
cd $FSFTUTDIR/my-study
cp ../fb1-raw-study/sensory-motor.par mgh-101.1/bold/003/sensory-motor.par
cp ../fb1-raw-study/sensory-motor.par mgh-101.1/bold/005/sensory-motor.par
cp ../fb1-raw-study/sensory-motor.par mgh-101.1/bold/007/sensory-motor.par
cp ../fb1-raw-study/sensory-motor.par mgh-101.1/bold/010/sensory-motor.par
}}}

While each of these have the same contents, they can just as easily be
different.

== Create a Session ID File ==

You would normally do the above steps for each of the subjects, but
you've probably got the point by now. If you were to create a session
for each, you would then have five sessions: mgh-101.1, mgh-103.1,
mgh-104.1, mgh-105.1, mgh-106.1. To help keep track of your sessions,
you can create a session ID file with your favorite text editor. This
is just a simple file with a list of your sessions. The contents should
look like (the order is unimportant):

{{{
mgh-101.1
mgh-103.1
mgh-104.1
mgh-105.1
mgh-106.1
}}}

This file can be passed to many of the FSFAST commands with the -sf
argument. Individual sessions can be passed with -s, as you will see
}}}

This is the Project Directory. You will run most of the FSFAST
commands from the Project Directory.

Run ls to see the contents of the Project Directory.
{{{
ls
}}}

You will see 18 folders with names like "sess09". These are the 18
subjects. There are some other files and folders there but don't worry
about them right now.

= Understanding the FS-FAST Directory Structure =

All of these sessions have been analyzed with the exception of
sess01.noproc. This session is what the directory structure
should look like immediately prior to beginning analysis. This
includes:

  * Directory stucture
  * Raw data
  * subjectname file
  * Parasigm files for each run

The directory structure and raw data are usually created by
"unpacking" the data with dcmunpack or unpacksdcmdir, but it could
also be done by hand. The subjectname file and paradigm files must be
added manually.

== The 'Session' Folder ==

The folder/directory where all the data for a session are stored is
called the 'session' or the 'sessid'. There may be more than one
session for a given subject (eg, in a longitudinal analysis). Go into
the sess01.noproc folder and run 'ls':

{{{
cd sess01.noproc
ls
}}}

You will see see two folders ('bold' and 'rest') and a file called
'subjectname'.

There is a file called 'sessidlist' that has a list of all the
sessions in this project. These types of files are called 'sessid
files'. They are convenient for

* Processing multiple sessions with one command-line
* Grouping sessions together
* Running analyses in parallel

== The 'subjectname' File ==

subjectname is a text file with the name of the FreeSurfer subject as
found in $SUBJECTS_DIR (ie, the location of the anatomical analysis).

View the contents of the text file (with 'cat', 'more', 'less',
'gedit', or 'emacs') and verify that this subject is in the
$SUBJECTS_DIR.

NOTE: it is important that the anatomical data and the functional data
be from the same subject. The contents of the subjectname file is
the only link! Make sure that it is right! There is a check for this
Line 359: Line 106:

---------------------------------------------------------------------
---------------------------------------------------------------------
= Preprocessing of fMRI Data (preproc-sess) =
== What preprocessing stages do you want to run? ==

There are at least seven types of preprocessing that are run on fMRI
data:

 * 1. motion correction (MC)
 * 2. spatial smoothing
 * 3. brain masking
 * 4. slice-timing correction (STC)
 * 5. intensity normalization
 * 6. resampling to common space
 * 7. B0 distortion correction.

Not all packages run all of these seven, and they are not always run
in the same order, and some stages are sometimes run as
post-processing. In FSFAST, we can run 1. MC, 2. smoothing 3. STC (optional),
and 4. intensity normalization (optional). We create a brain mask,
but we do not mask the functional data.

For these exercises, cd into the fb1-raw-study directory:
{{{
cd $FSFTUTDIR/my-study
}}}
Note that the fully preprocessed data can be found in fb1-preproc-raw.

== Run preprocessing ==

To run MC and spatial smoothing by 5 mm FWHM along with brain mask
creation on one session run:

{{{
cd $FSFTUTDIR/my-study
preproc-sess -s mgh-101.1 -fwhm 5
}}}

Note that this has already been done with all subjects in the
fb1-preproc-study directory. To preprocess all of the session, you
would run "preproc-sess -sf sessid -fwhm 5". Also note that if
preprocessing has already been performed on a session, it will
automatically skip it and move on to the next (unless you add -force
to the command-line).

== Examine additions to FSFAST hierarchy ==

As mentioned above, each processing step adds something to the FSFAST
hierarchy. To see what has been created, run:

{{{
cd $FSFTUTDIR/my-study
ls mgh-101.1/bold/003
}}}

You will see:
 * f.nii (the raw data)
 * fmc.nii (motion corrected)
 * fmcsm5.nii (motion corrected and smoothed).
 * fmc.mcdat - text file with the motion correction parameters (AFNI)

You will also see mcextreg.bhdr. This is a binary file with the
orthogonalized motion correction parameters which can later be used as
nuisance regressors when analyzing the data. These files will exist in
each of the runs (i.e., 005, 007, 010).

{{{
cd $FSFTUTDIR/my-study
ls mgh-101.1/bold/masks
}}}

You will see a brain.nii volume in the masks folder. This is a binary
mask of the brain as found by FSL's BET program. The functional data
themselves are not masked.<<BR>>

To view the translation components, run
{{{
cd $FSFTUTDIR/my-study
plot-twf-sess -s mgh-101.1 -mc
}}}
This will bring up a plot with the translations for each of the
runs. Notice that run 007 has a head jerk at about time point 11.

---------------------------------------------------------------------
---------------------------------------------------------------------
= Function-Structure Registration =

In order to render the functional results on the anatomical background
as well as to map the functional results into a common space for group
analysis, it is necessary to register/align the functional volume with
a structural volume. In FSFAST, we:

 * First register to the same-subject FreeSurfer anatomical with a 6 DOF registration.
 * Map the functional to Talairach/MNI305/fsaverage space by concatenating the within-subject function/structure registration with the Talairach registration (talairach.xfm) created when the subject was reconstructed.
 * For mapping to surface-based space, we concatenate the within-subject function/structure registration with the surface-based registration.

Since we are only dealing with the functional analysis here, we will
just consider the within-subject function/structure registration.

For these series of exercises, cd into fb1-preproc-study:
{{{
cd $FSFTUTDIR/fb1-preproc-study
}}}

== View unregistered (tkregister-sess) ==

Run the following command to see how the functional and structural are
aligned prior to performing any automatic registration (delete the old one, if there).

{{{
cd $FSFTUTDIR/fb1-preproc-study
rm -f mgh-101.1/bold/register.dat
tkregister-sess -s mgh-101.1 -regheader -inorm
}}}

You should see the anatomical image (left, below). The green line is
the cortical surface for that subject. Hit the compare button to see
the functional middle image. The cortical surface shows that the
volumes are not in line.

{{attachment:tkregister-anat.gif}}
{{attachment:tkregister-before.gif}}
{{attachment:tkregister-after.gif}}

== Run automatic registration (spmregister-sess) ==
{{{
cd $FSFTUTDIR/fb1-preproc-study
fslregister-sess -s mgh-101.1
}}}
When this command is complete, you will see a register.dat file in
mgh-101.1/bold. This is the only change. The functional data are not
resampled! Note that this command requires FSL. You can also run:
{{{
cd $FSFTUTDIR/fb1-preproc-study
spmregister-sess -s mgh-101.1
}}}
This does the same thing but requires matlab and spm in
your matlab path (in startup.m), otherwise use fslregister-sess.

== Check automatic registration (tkregister-sess) ==
{{{
cd $FSFTUTDIR/fb1-preproc-study
tkregister-sess -s mgh-101.1 -inorm
}}}
When you hit the Compare button, you should see the image on the
right, above.

== Check talairach registration ==

When performing a volume-based group analysis, it is important that
the talairach transform be accurate (or as accurate as talairach can
be). You should have done this during the FreeSurfer reconstruction
for this subject. The command below is just a reminder of how to do
it.

{{{
tkregister2 --s fbph1-101 --fstal --surf
}}}

---------------------------------------------------------------------
---------------------------------------------------------------------
= First-Level Analysis =

The First-Level Analysis (FLA) consists of:

 * Setting up models of the task-related components
 * Setting up models of the nuisance components
 * Defining contrasts
 * Fitting the model
 * Making inferences.

In FSFAST, the FLA is done in two stages. In the first stage, the
FLA is configured (with mkanalysis-sess). This is done ONCE regardless
of how many data sets you have (you do not even need to have any data
to run the configuration).

In the second stage, you actually perform the analysis/inference with
selxavg3-sess by passing it the configuration and the session that you
want to analyze. selxavg3-sess customizes the analysis for that
session based on what it finds in the hierarchy, builds the design
matrix and performs the analysis. Breaking the FLA up into these two
stages assures that all sessions are analyzed in the same way.

In this exercise, we will configure two analyses, one assuming an HRF
and the other using a Finite Impulse Response (FIR) model. You will
then launch the analysis and view the results. The fully analyzed data
can be found in fb1-analysis-study. These exercises will be done in
fb1-preproc-study:

{{{
cd $FSFTUTDIR/fb1-preproc-study
}}}

== Configure Analysis and Contrasts I: Gamma HRF Model ==

Configuring the FLA is performed with mkanalysis-sess. When you run:
{{{
cd $FSFTUTDIR/fb1-preproc-study
mkanalysis-sess -gui
}}}

You will see the following window:

{{attachment:mkana.init.gif}}

You will use this window to specify the input of the analysis, the
hemodynamic response model, contrasts and nuisance regressors.
The red fields are fields that you must enter before you can save
the analysis. There is a lot going on with this GUI, so we'll break it
down. Note that many of the components have "tooltips" that will show
when you pause the mouse pointer over them.

In the upper left corner is a panel called "FS-FAST Hierarchy". The
"Func Stem" is the input to the analysis. You should specify the
output from the preprocessing. For this exercise, we are going to use
the motion corrected and 5mm smoothed data. This functional volume is
called fmcsm5.nii in the hierarchy which makes its stem "fmcsm5" (i.e.,
just strip off the nii). Enter "fmcsm5" into the field. When you hit
return, it changes from red to white. Next, enter the TR (sec). For
this experiment it was 3 sec. This will be checked against the TR
found in the input nifti file. Leave INorm checked.
{{{
Func Stem = fmcsm5
TR = 3
}}}

Turn your attention to the "Noise and Nuisance Variables" panel. Low
frequency noise so prevalent in fMRI is compensated for in a
combination of three ways. Drift components are modeled with
polynomial regressors. The order can be adjusted, but leave it at 2
for now. The motion correction parameters can be used as regressors by
checking the "MC Regressors" box (do so now). Finally, the remaining noise is
modeled as time-invariant linear AR1 process when the "Temporal
Whitening" box is checked (leave it so). There is one additionally way
to compensate for noise through the use of a "Time Point Exclude
File", but we will not consider that here.
{{{
Check "MC Regressors" Box
}}}

You will specify the model of the task-related signal in the "Event
Related/Block Design" panel (leave that box checked). Choose the
number of conditions by clicking on the "NConditions" slider. This is
the number of TASK conditions (do not include the Null/Fixation
condition). In this example, we have two conditions (Odd and Even), so
set this to 2. To the right of this is the "Paradigm File". Enter
"sensory-motor.par". Note that the number of task conditions in the
paradigm file must match that specified with "NConditions". Below, you
will specify the Hemodynamic Response Model. There are three choices:
Gamma, SPM HRF, and FIR. Choose Gamma for now. If you hit the "Plot"
button it will show the Gamma and SPM HRF. As you change the Gamma
parameters (Delay, Dispersion (Tau), and Exponent (Alpha)), the Gamma
plot will change. Make sure that they are at Delay=2.25, Tau=1.25, and
Alpha=2.
{{{
Set NConditions to 2
Paradigm File = sensory-motor.par
Press "Plot" button
}}}

At this point, you have specified the model of the BOLD signal
including HRF, nuisance, and noise. The GUI should look like the image
below:

{{attachment:mkana-gamma-precon.gif}}

Now you are ready to specify contrasts. A contrast is an instantiation
of a hypothesis and is represented by a contrast matrix (i.e., a linear
summation of the regression coefficients). Contrasts are managed
through a separate GUI accessed through the "Contrast" list box. When
you click on "Add Contrast", you will see the following screen:
{{{
Click on "Add Contrast"
}}}

{{attachment:mkcon.init.gif}}

There are several things going on here, but the most important is the
list of conditions in the middle of the GUI (i.e., "Condition 1",
"Condition 2") shown as green, red, and black radio buttons.

 * Green indicates an "active" condition.
 * Red means a "control" condition
 * Black means to ignore the condition in the contrast.

Active conditions are given a weight of +1; controls are given -1;
ignores get 0. The weight is given to the right of the buttons. All
contrasts are implicitly computed against the Null or Fixation
condition. If you want to test the null hypothesis that Condition 1 is
no different than the Null condition, then you would make Condition 1
active and ignore the rest. To test the null hypothesis that Condition
1 is no different than Condition 2, then you would make Condition 1
active and Condition 2 control.

For this exercise, we are going to test four NULL hypotheses:

 * Odd == Fixation (odd-v-0)
 * Even == Fixation (even-v-0)
 * Odd == Even (odd-v-even)
 * Odd+Even == Fixation (odd+even)

The last one tests whether the average of the responses to odd and
even are different than fixation. Remember that, according to the
Paradigm File, Condition 1 is Odd, and Condition 2 is Even.

{{{
Click on the green button next to "Condition 1". (see its weight change from 0 to 1).
}}}

Condition 2 is already black (ignored). This our first contrast, so there is
nothing more we need to do except give the contrast a name. You should give
your contrasts meaningful but terse names.
{{{
Specify "odd-v-0" for the name of this contrast.
Hit the "Done/Save" button.
}}}
You will now see "odd-v-0" appear in the Contrast list box in the mkanalyiss GUI.

{{{
Click on "Add Contrast" again to bring up the contrast GUI again.
Click on the green button next to Condition 2 (see its weight change from 0 to 1).
}}}

This activates condition 2. Condition 1 is still black (ignored).
{{{
Change the name to "even-v-0"
Click Done/Save
}}}
"even-v-0" will appear in the list box.

{{{
Click on "Add Contrast" again
Click the green button next to Condition 1 (weight changes to 1.0).
Click the red button next to Condition 2 (weight changes to -1.0).
Change the name to "odd-v-even"
Click Done/Save.
}}}

{{{
Click on "Add Contrast" one more time
Click the green button next to Condition 1. Weight changes to 1.0.
Click the green button next to Condition 2. Both weights change to 0.5.
Change the name to "odd+even"
Click Done/Save.
}}}

You can go back and view and/or edit an contrast by clicking on it in
the list box.

The last thing you have to do is to give your analysis a name. Like
the contrasts, it should be terse but descriptive. It cannot have any
spaces or blanks.
{{{
Specify "sm-gamma-fwhm5" as the analysis name.
}}}
sm = sensory-motor, gamma = Gamma HRF, and fwhm5 for the input. The interface should now look

{{attachment:mkana.final.gif}}

{{{
Hit the "Save" button
Hit the "Quit" button
}}}

After you hit Quit, control will be returned to the shell that you ran
mkanalysis-sess from. If you type "ls", you will see a new folder
called "sm-gamma-fwhm5". If you "ls sm-gamma-fwhm5", you will see
analysis.info, analysis.cfg, odd-v-0.mat, even-v-0.mat,
odd-v-even.mat, odd+even.mat. Your configuration is stored in these
files. You can browse/edit your configuration by running:

{{{
mkanalysis-sess -gui -analysis sm-gamma-fwhm5
}}}

== Configure Analysis and Contrasts II: FIR HRF Model ==

Now we are going to use a Finit Impulse Response (FIR) to model the
hemodynamic response. The FIR does not make any assumptions about the
shape of the HRF but is also less interpretable. Again, run

{{{
mkanalysis-sess -gui
}}}

Set the Func Stem, TR, NConditions, and Paradigm File as above, but
now click on the "FIR" checkbox. This will enable the "Total Time
Window", "PreStim", and "TER" entry boxes. The Time Window is the
window within which we will estimate the HRF. Given that the task is
15 sec long and the rest is 15 sec, let's choose 27 sec. The PreStim
is the amount of time before stimulus onset to start estimating the
HRF. A non-zero PreStim gives us an idea of what the baseline is at
stimulus onset. Set it to 6.

Set up the same contrasts as you did above, then name the analysis
"sm-fir-fwhm5", hit Save, then Quit.

== Analyze First Level (selxavg3-sess) ==

You are now ready to analyze some data! Note that the fully analyzed
data (along with correctly configured analyses) can be found in
fb1-analysis-study. To analyze the data for session mgh-101.1 with the
sm-gamma-fwhm5 analysis, run:

{{{
selxavg3-sess -s mgh-101.1 -analysis sm-gamma-fwhm5
}}}

Note that if you want to analyze all the sessions, you can run
"selxavg3-sess -sf sessid -analysis sm-gamma-fwhm5", but this might
take a while.

To analyze the data for session mgh-101.1 with the sm-fir-fwhm5
analysis, run:
{{{
selxavg3-sess -s mgh-101.1 -analysis sm-fir-fwhm5
}}}

== Examine additions to the FSFAST hierarchy ==

{{{
ls mgh-101.1/bold
ls mgh-101.1/bold/sm-gamma-fwhm5
}}}

The bold folder now has sm-gamma-fwhm5 (as well as the fir too). This
folder has all of the results for this analysis, including:

 * beta.nii - regression coefficients
 * rvar.nii - residual error variance
 * mask.nii - mask (copy of bold/masks/brain.nii)
 * meanfunc.nii - mean functional image
 * fsnr.nii - functional SNR map
 * X.mat - design matrix (in matlab format)
 * dof - text file with degrees of freedom
 * fwhm.dat - text file with smoothness estimate (Full-Width/Half-Max)
 * contrast folders: odd-v-0, even-v-0, odd+even, odd-v-even

Look in one of the contrasts:

{{{
ls mgh-101.1/bold/sm-gamma-fwhm5/odd-v-0
}}}

 * ces.nii - contrast effect size (contrast matrix * regression coef)
 * cesvar.nii - variance of CES
 * sig.nii - significance map (-log10(p))

== Visualize ==

Visualization of the first-level analysis will be done in the folder
with the fully analyzed data set:

{{{
cd $FSFTUTDIR/fb1-analysis-study
}}}

=== Volume-based visualization (tkmedit-sess) ===

View the result of the Gamma HRF analysis on the FreeSurfer anatomical
volume for mgh-101.1 with:
{{{
tkmedit-sess -s mgh-101.1 -aparc+aseg -analysis sm-gamma-fwhm5 \
  -c odd-v-0 -c even-v-0 -c odd+even -c odd-v-even
}}}
This will automatically load sig.nii.
When you configure the functional overlay with
View->Configure->Functional Overlay, you will see that there are 4
"Time Points". Each point is different contrast (i.e., 0 is odd-v-0, 1
is even-v-0, etc). Scroll through each one.

View the result of the HRF HRF analysis for mgh-101.1 with:
{{{
tkmedit-sess -s mgh-101.1 -aparc+aseg -analysis sm-fir-fwhm5 -c odd-v-0
}}}

Here we view only one contrast at a time because each contrast has
multiple time points for each point in the FIR time window. When you
click on a point, you will see the HRF for both Condition 1 (Odd) and
Condition 2 (Even) blocks.

Finally, view the overlay maps from the Gamma with the HRF from the FIR:
{{{
tkmedit-sess -s mgh-101.1 -aparc+aseg -analysis sm-fir-fwhm5 \
  -mapanalysis sm-gamma-fwhm5 -c odd-v-0 -c even-v-0 -c odd+even -c odd-v-even
}}}
When you configure the overlay, you will see that there are 4 "Time
Points" -- these correspond to the 4 contrasts from the Gamma analysis.

=== Surface-based visualization (tksurfer-sess) ===
View the result of the Gamma HRF analysis on the FreeSurfer anatomical
surface for mgh-101.1 with:
{{{
tksurfer-sess -hemi lh -aparc -s mgh-101.1 -analysis sm-gamma-fwhm5 \
  -c odd-v-0 -c even-v-0 -c odd+even -c odd-v-even
}}}
Note that the command-line is nearly identical to that of tkmedit
above. The difference is that the hemisphere is specified ("-hemi
lh"), and "-aparc+aseg" is replaced with "-aparc" to load the
surface-based segmentation.

---------------------------------------------------------------------
---------------------------------------------------------------------

= Higher-Level (Group) Analysis =

Higher-Level is where you make inferences about the population that
your subjects are drawn from. It is a bit confusing at times because
both use GLMs, so at both levels you are constructing design matrices,
contrasts, etc.

Traditionally, fMRI group analysis has been done in a standard volume
space (i.e., Talairach/MNI152/MNI305). With FreeSurfer, we also have the
option to analyze group data in the surface space. Volume-based
analyses are done in MNI305 space (which is the same as the fsaverage
subject).

For the next exercises, we will work in the fb1-analysis-study directory
{{{
cd $FSFTUTDIR/fb1-analysis-study
}}}
All the first level analyses have been done here. The group-level
analyses have been done in the group-analysis-tut directory.

== Volume-based (MNI305/fsaverage) Group Analysis ==

=== Assemble the Data (isxconcat-sess) (Volume) ===

The first step in the group analysis is to "assemble" the data. This
means creating a single 4D file with where the 4th "time" dimension is
actual all the subjects concatenated together in a common space. There is a
different command, depending upon whether the common space is volume-
or surface-based.

To run the volume-based concatenation, run the command below. Note
that the data from this command already exist in the
group-analysis-tut directory.

{{{
isxconcat-sess -sf sessid -analysis sm-gamma-fwhm5 -c odd-v-0 -o group-analysis
}}}

This command will:
 * Go through each session in the sessid file
 * Find the odd-v-0 contrast in the sm-gamma-fwhm5 analysis
 * Use the register.dat for that session to resample to MNI305/fsaverage space.
 * Concatenate all subjects' data into one file

The concatenated data are saved in
group-analysis/sm-gamma-fwhm5/odd-v-0/tal.ces.nii file. You can run
other contrasts by adding -c arguments. In addition to this output,
several other files are created. To see them,

{{{
cd $FSFTUTDIR/fb1-analysis-study
ls group-analysis-tut/sm-gamma-fwhm5
cat group-analysis-tut/sm-gamma-fwhm5/sessid.txt
cat group-analysis-tut/sm-gamma-fwhm5/ffxdof.dat
}}}

In the output directory, you will see a series of files that start
with "tal":
 * tal.meanfunc.nii is a stack where each "time point" is the mean functional image of each subject sampled in the MNI305 space.
 * tal.masks.nii are the binary masks for all the subjects
 * tal.fsnr.nii are the functional SNR maps from each subject.
 * tal.mask.nii is a single binary mask made from the intersection of the individuals.
 * ffxdof.dat is the fixed-effects DOF across all subjects.
 * sessid.txt is the list of sessions, the corresponding freesurfer subject name, and the DOF contributed by each subject.

You will also see some files that being with "lh". These are the same
thing in the surface-based space (see below).

Now look in the directory for odd-v-0 the contrast
{{{
cd $FSFTUTDIR/fb1-analysis-study
ls group-analysis-tut/sm-gamma-fwhm5/odd-v-0
}}}

You will see:
 * tal.ces.nii - the contrast maps for each of the subjects
 * tal.cesvar.nii are the variance of the contrast for each subject (i.e., the square of the standard error). This variance is needed for fixed-effects and weighted random-effects analysis.


=== Quality Assurance (Volume) ===

There are three important quality assurance steps that can be performed
here. First, view the mean functionals to make sure that all are
registered together properly. To do this run,

{{{
tkregister2 --s fsaverage --surf --regheader --check-reg \
  --mov group-analysis-tut/sm-gamma-fwhm5/tal.meanfunc.nii
}}}
The image window will show the MNI305 brain. Hit the "Compare" button
to show the average functional of the first session. Click in the
image window, then hit the 'a' key. Each time you hit the 'a' key, it
will advance to the next subject.

The next QA step is to check the individual masks. This can be done
     == Functional Subdirectories (FSDs) ==

The other two directories (bold and rest) are 'functional
subdirectories' (FSDs) and contain functional data. If you 'ls rest'
you will see '001'. If you 'ls bold' you will see '001 002 003
004'. Each of these is a 'run' of fMRI data, ie, all the data
collected from a start and stop of the scanner.

=== Raw Data ===

Go into the first run of the bold directory with 'cd bold/001' and run
ls. You will see 'f.nii.gz wmfir.par workmem.par'. The raw data is
stored in f.nii.gz (compressed NIFTI) and is directly converted from
the DICOM file. Examine this file with mri_info:

{{{
mri_info --dim f.nii.gz
mri_info --res f.nii.gz
}}}

The first command results in '64 64 30 142'. This is the dimension of
the functional data. Since it is functional, it has 4 dimensions: 3
spatial and one temporal (ie, 64 rows, 64 cols, 30 slices, and 142
time points or TRs or frames). The second command results in '3.438
3.437 5.000 2000.000'. This is the resolution of the data, ie, each
voxel is 3.438mm by 3.437mm by 5.000mm and the TR is 2000ms (2 sec).

View the functional data with:
{{{
tkmedit -f f.nii.gz -t f.nii.gz
}}}
Click on a point to view the waveform at that point.


=== Paradigm Files ===

The workmem.par and wmfir.par files are paradigm files. They are text
files that you create that indicate the stimulus schedule (ie, which
stimulus was presented when).

View workmem.par. Each row indicates a stimulus presentation. You will
see that each row has 5 columns. The columns are:
  
  * Stimulus Onset Time (sec)
  * Numeric Stimulus Identifier
  * Stimulus Duration (sec)
  * Weight (usually 1)
  * Text Stimulus Identifier (redundant with Numeric Stimulus Identifier)

The Stimulus Onset Time is the onset relative to the acquiistion time
of first time point in f.nii.gz. The Numeric and Text Stimulus
Identifiers indicate which stimulus was presented. The Stimulus
Duration is the amount of time the stimulus was presented.

In this case, there are 5 event types:

  * Encode - encode phase
  * EDistrator - emotional distractor
  * EProbe - probe after emotional distractor
  * NDistrator - neutral distractor
  * NProbe - probe after neutral distractor

Note two things: (1) Not all the time is taken up, and (2)
Basline/Fixation is not explicitly represented. By default, any time
not covered by stimuluation is assumed to be baseline.

  * What time was the third Encode presented?
  * What time was the last Neutral Distractor presented?
  * Is this timing for all runs the same?


= Preprocessing =

Once the data have been arranged in the above directory structure and
naming convention, they are ready to be preprocessed. In FS-FAST, it
is assumed that each data set will be analyzed in three ways:

  * Left Cortical Hemisphere
  * Right Cortical Hemisphere
  * Subcortical Structures

You will need to decide how much to smooth the data and whether you
want to do slice-timing correction. In this analysis, we will smooth
the data by 5mm Full-Width/Half-Max (FWHM) and correct for slice
timing. The slice-timing for this particular data set was 'Ascending',
meaning that the first slice was acquired first, the second slice was
acquired second, etc. To preprocess the data, run:

{{{
preproc-sess -s sess01 -surface fsaverage lhrh -mni305 \
   -fwhm 5 -stc up -fsd bold -per-run
}}}

This data has already been preprocessed, so it should just verify that
it is up-to-date and return. This command has several arguments:

  * -s sess01 : indicates which session to preprocess
  * -surface fsaverage lhrh : indicates that data should be sampled to the left and right hemispheres of the 'fsaverage' subject
  * -mni305 : indicates that data should be sampled to the mni305 volume (2mm isotropic voxel size)
  * -fwhm 5 : smooth by 5mm FWHM
  * -stc up : slice-timing correction with ascending ('up') slice order
  * -fsd bold : indicate the functional subdirectory
  * -per-run : register each run separately

This command does a lot (and it can take quite a long time to run). To
understand what it does, we will go back into one of the run
directories and see what it creates. To do this, 'cd sess01/bold/001'
and type 'ls'. This directory previously held only f.nii.gz,
workmem.par, and wmfir.par; now there are a lot of files, each
indicative of a different preprocessing stage. Now type 'ls -ltr'. This
command sorts the files by creation time with the oldest at the top
and the newest at the bottom. The preprocessing is progressive,
meaning that the output of one stage is the input to the next.

== Template ==

This stage creates template.nii.gz (and template.log). This is the
middle time point from the raw functional data (f.nii.gz). This is
the reference used to motion correct and register the functionals
for this run to the anatomical. It is also used to create masks of the
brain.

* Verify that it has the same dimension and resolution as the raw data using mri_info.

== Masking ==

The masks for this run are stored in the 'masks' directory. Run 'ls
-ltr masks'. You will see a file called 'brain.nii.gz'. This is a
binary mask created using the FSL BET program. There is also a file
called 'brain.e3.nii.gz' which is the mask eroded by three
voxels. These have the same dimensions as the template. View the masks
Line 960: Line 239:
{{{
tkmedit fsaverage orig.mgz \
  -overlay group-analysis-tut/sm-gamma-fwhm5/tal.masks.nii -fthresh 0.5
}}}
The threshold of 0.5 is appropriate because these masks are binary
(ie 0-1). When you View->Configure->Functional Overlay, you will see
that there are 5 "Time Points" (0-4) corresponding to the 5 subjects.
Advance through each one to assure that the masks are in the proper place.

The final step is to look at the functional SNR maps with
{{{
tkmedit fsaverage orig.mgz \
  -overlay group-analysis-tut/sm-gamma-fwhm5/tal.fsnr.nii \
  -timecourse group-analysis-tut/sm-gamma-fwhm5/tal.fsnr.nii \
  -fthresh 50 -fmax 250
}}}

{{{
When you View->Configure->Functional Overlay, you will see
that there are 5 "Time Points" (0-4) corresponding to the 5 subjects.
Advance through each one to assure that the SNR maps are "consistent".
}}}

When you click on a voxel, you will see the SNR for each subject
plotted in the "Time Course" window. The actual value of the FSNR will
vary depending upon how much smoothing you've done and the details of
the acquisition. You are looking for outliers here.

=== Group GLM Analysis (Volume) ===

When you perform a group analysis, you are looking for effects of the
task across the population. In this example, we are only going to test
whether the mean effect of Odd vs Fixation (odd-v-0) is different than
0. This is known as a one-sample-group-mean (OSGM) and corresponds to
a group design matrix that is simply a column of 1s. One can perform
considerably more elaborate tests (e.g., differences in groups). Now, cd
to where the data are:

{{{
cd $FSFTUTDIR/fb1-analysis-study/group-analysis-tut/sm-gamma-fwhm5/odd-v-0
}}}

==== Random Effects (RFx, OLS) (Volume) ====
Random effects is a test of whether the mean of the population that
the subjects were drawn from is 0. It is done with mri_glmfit:

{{{
cd $FSFTUTDIR/fb1-analysis-study/group-analysis-tut/sm-gamma-fwhm5/odd-v-0
mri_glmfit --y tal.ces.nii --osgm --mask ../tal.mask.nii --glmdir tal.rfx.osgm --nii
}}}

 * tal.ces.nii are the inputs for this contrast for all subjects.
 * --osgm tells it to create the simple OSGM design matrix and to create a group contrast to test the group mean against 0.
 * The mask is used to exclude extra-brain voxels.

The results will be stored in tal.rfx.osgm:
{{{
cd $FSFTUTDIR/fb1-analysis-study/group-analysis-tut/sm-gamma-fwhm5/odd-v-0
ls tal.rfx.osgm
ls tal.rfx.osgm/osgm
}}}

You will see several files in the output. beta.nii is the map of
regression coefficients. There is also an "osgm" folder with
sig.nii. This is a map of the significances of the test. We will view
it below.

==== Weighted Random Effects (WRFx, WLS) (Volume) ====
In a real experiment, some subjects are noisier than others, and it is
a good idea to take this into account since we have information about
how noisy a subject is through the lower-level analysis. In weighted
least squares (WLS), this is handled by weighting each subject by the
inverse of their noise (i.e., noisier subjects get lower weight). To do
this with mri_glmfit, run:

{{{
cd $FSFTUTDIR/fb1-analysis-study/group-analysis-tut/sm-gamma-fwhm5/odd-v-0
mri_glmfit --y tal.ces.nii --osgm --glmdir tal.wrfx.osgm --nii --mask ../tal.mask.nii \
    --wls tal.cesvar.nii
}}}

The command-line is almost the same as the RFx except that the
first-level noise variances (tal.cesvar.nii) are passed with the --wls
option. The results are saved in tal.wrfx.osgm and have the same
naming convention as the RFx.

==== Fixed Effects (FFx) (Volume) ====
A fixed effects analysis tests whether an effect exists within the
given subject pool (as opposed to the population that the subjects
were pulled from). It is like doing one big fist-level analysis. The
big difference with RFx or WRFx is that the variance is computed from
the lower level analysis variance. This is done in mri_glmfit with:

{{{
cd $FSFTUTDIR/fb1-analysis-study/group-analysis-tut/sm-gamma-fwhm5/odd-v-0
mri_glmfit --y tal.ces.nii --osgm --glmdir tal.ffx.osgm --nii --mask ../tal.mask.nii \
   --yffxvar tal.cesvar.nii --ffxdofdat ../ffxdof.dat
}}}

Again, the command-line is similar to above, except that the lower
level variances are passed with --yffxvar, and the total lower level
DOFs are passed with --ffxdofdat.

=== Output and visualization (Volume) ===

To visualize these three analyses, we will first concatenate them into
one file to make it easier to view in tkmedit:

{{{
cd $FSFTUTDIR/fb1-analysis-study/group-analysis-tut/sm-gamma-fwhm5/odd-v-0
mri_concat tal.rfx.osgm/osgm/sig.nii \
           tal.wrfx.osgm/osgm/sig.nii \
           tal.ffx.osgm/osgm/sig.nii \
        --o all.sig.nii
tkmedit fsaverage orig.mgz -aux brain.mgz \
     -overlay all.sig.nii -fthresh 2 -fmax 4
}}}
Configure the overlay with View->Configure->Functional Overlay. Scroll
through the "time" points to see the differences in the analyses.

=== Correction for Multiple Comparisons/Cluster Analysis (Volume) ===

With so many voxels in fMRI maps, it is very likely that many voxels
will appear to be active purely by random chance (ie, a false
positive). The is known as the "Problem of Multiple Comparisons". One
way around this is to do a cluster analysis in which active voxels are
eliminated unless they appear in a cluster, the idea being that false
positives will not appear next to each other. To do this in FSFAST,
run the following command on the RFx analysis:

{{{
cd $FSFTUTDIR/fb1-analysis-study/group-analysis-tut/sm-gamma-fwhm5/odd-v-0
mri_volcluster --in tal.rfx.osgm/osgm/sig.nii --thmin 2 \
   --fwhmdat tal.rfx.osgm/fwhm.dat --cwpvalthresh 0.1 \
   --mask tal.rfx.osgm/mask.nii --fsaverage \
   --sum tal.rfx.osgm/osgm/cluster.sum --sign pos\
   --cwsig tal.rfx.osgm/osgm/cwsig.cluster.nii \
   --out tal.rfx.osgm/osgm/sig.cluster.nii \
   --ocn tal.rfx.osgm/osgm/ocn.cluster.nii
}}}

In this command we use the significance map as input thresholded at 2,
where 2 is -log10(p), so p < .01. The --fwhmdat option tells
mri_volcluster to use the spatial smoothness estimate from GLM
analysis (saved in tal.rfx.osgm/fwhm.dat). The "--cwpvalthresh 0.1"
tells it to ignore clusters unless their p-value is less than 0.1 (you
can threshold more stringently later on).

One of the outputs is a cluster table (cluster.sum). This is a text
file with the list of clusters. View it with:

{{{
tkmedit -f template.nii.gz -overlay masks/brain.nii.gz -fthresh 0.5
tkmedit -f template.nii.gz -overlay masks/brain.e3.nii.gz -fthresh 0.5
}}}

The brain.nii.gz is used to constrain voxel-wise operations. The
eroded mask (brain.e3.nii.gz) is used to compute the mean functional
value used for intensity normalization and global mean time
course. There are other masks there that we will get to later.

== Intensity Normalization and Global Mean Time Course ==

By default, FSFAST will scale the intensities of all voxels and time
points to help assure that they are of the same value across runs,
sessions, and subjects. It does this by dividing by the mean across
all voxels and time points inside the brain.e3.nii.gz mask, then
multiplying by 100. This value is stored in global.meanval.dat. This
is a simple text file which you can view. At this point, this value is
stored and used later. A waveform is also constructed of the mean at
each time point (text file global.waveform.dat). This can be used as a
nuisance regressor.

* What are the global means for runs 1, 2, 3, and 4?

== Functional-Anatomical Cross-modal Registration ==

The next six files (init.register.dof6.dat, register.dof6.dat,
register.dof6.dat.mincost, register.dof6.dat.sum,
register.dof6.dat.log, register.dof6.dat.param) deal with the
registration from the functional to the same-subject FreeSurfer
anatomical. There are only two files here that are really
important: register.dof6.dat and register.dof6.dat.mincost.

 * register.dof6.dat - is a text file that contains the registration matrix.
 * register.dof6.mincost - is a text file that contains a measure of
     the quality of the registration.

The registration is will be revisited below when we talk about Quality
Assurance

== Motion Correction (MC) ==

The motion correction stage produces these files: fmcpr.mat.aff12.1D,
fmcpr.nii.gz, mcprextreg, mcdat2extreg.log, fmcpr.nii.gz.mclog,
fmcpr.mcdat. There are only three important file here:

  * fmcpr.nii.gz -- this is the motion corrected functional data. It has the same size and dimension as the raw functional data.
  * fmcpr.mcdat - text file of the amount of motion at each time point. This is important for Quality Assurance (see below).
  * mcprextreg - text file of the motion correction parameters assembled into an orthogonalized matrix that can be used as nuisance regressors

  * Verify that fmcpr.nii.gz has the same dimension as f.nii.gz using mri_info.

== Slice-Timing Correction (STC) ==

Slice-timing corretion compensates for the fact that each of the 30
slices was acquired separately over the course of 2 sec. It does this
by interpolating between time points to align each slice to the time
of the middle of the TR. The file created with this is fmcpr.up.nii.gz
(and fmcpr.up.nii.gz.log).

  * Verify that fmcpr.up.nii.gz has the same dimension as f.nii.gz using mri_info.

== Resampling to Common Spaces and Spatial Smoothing ==

At this point, the functional data has stayed in the 'native
functional space', ie, 64x64x30, 3.4x3.4x5mm3. Now it will be sampled
into the 'Common Space'. The Common Space is a geometry where all
subjects are in voxel-for-voxel registration. There are three such
spaces in FSFAST:

  * Left hemisphere of fsaverage (fmcpr.up.sm5.fsaverage.lh.nii.gz)
  * Right hemisphere of fsaverage (fmcpr.up.sm5.fsaverage.rh.nii.gz)
  * Volume of fsaverage (MNI305 space) - for subcortical analyses (fmcpr.sm5.mni305.2mm.nii.gz)

Each of these is the entire 4D functional data set resampled into the
common space. The spatial smoothing is performed after
resampling. Surface-based (2D) smoothing is used for the surfaces; 3D
for the volumes.

Check the dimensions of the MNI305 space volume:
{{{
mri_info --dim fmcpr.up.sm5.mni305.2mm.nii.gz
mri_info --res fmcpr.up.sm5.mni305.2mm.nii.gz
}}}
The dimension will be '76 76 93 142' meaning that there are 76
columns, 76 rows, 93 slices but still 142 time points (same as the raw
data). The resolution will be '2.0 2.0 2.0 2000' meaning that each
voxel is 2mm in size and the TR is still 2 sec. The transformation
to this space is based on the 12DOF talairach.xfm created during the
FreeSurfer reconstruction.

Check the dimensions of the Left Hemisphere 'volume':
{{{
mri_info --dim fmcpr.up.sm5.fsaverage.lh.nii.gz
mri_info --res fmcpr.up.sm5.fsaverage.lh.nii.gz
}}}

The dimension is '163842 1 1 142'. This 'volume' has 163842 'columns',
1 'row', and 1 'slice' (still 142 time points). You are probably
confused right now. That's ok, it's natural. At this point the notion
of a 'volume' has been lost. Each 'voxel' is actually a vertex (of
which there are 163842 in the left hemisphere of fsaverage). Storing
it in a NIFTI 'volume' is just a convenience.

The 'resolution' is '1.0 1.0 1.0 2000'. The values for the first 3
dimensions are meaningless because there are no columns, rows, or
slices on the surface so the distances between them are
meaningless. The last value indicates the time between frames and is
still accurate (2 sec).

The transformation to this space is based on the surface-based
intersubject registration created during the FreeSurfer
reconstruction.

== Quality Assurance ==

=== Motion Correction ===

The motion plots can be viewed with:
Line 1112: Line 361:
more tal.rfx.osgm/osgm/cluster.sum plot-twf-sess -s sess01 -fsd bold -mc
Line 1115: Line 364:
You can also see it here: FsFastTutorial/ClusterSummary. The table
shows the size of each cluster in voxels (mm^3), the Talairach
coordinate, the maximum significance in the cluster, and GRFCWP (the
clusterwise p-value computed using gaussian random field). Some of them are 0, meaning that the
p-value was so low that it could not be computed (remember, this is a
good thing!).

The other outputs are volumes. cwsig.cluster.nii is a map of the
clusters with the voxel value equal to the -log10(pvalue) of the
cluster that the voxel is associated with. sig.cluster.nii is the
original sig map with non-cluster voxels removed. ocn.cluster.nii is a
map where the value at each voxel is replaced by the number of the
cluster that the voxel is associated with. View the clusterwise
significance:

{{{
tkmedit fsaverage orig.mgz -aux brain.mgz \
     -overlay tal.rfx.osgm/osgm/cwsig.cluster.nii -fthresh 2 -fmax 4
}}}

--------------------------------------------------------------------
--------------------------------------------------------------------

== Surface-based Group Analysis ==

=== Assemble the Data (isxconcat-sess) (Surface) ===

The surface-based analysis proceeds in much the same way as the
volume-based. The assembly command is the same except you add "-hemi
lh" (or "-hemi rh") to the command:
{{{
cd $FSFTUTDIR/fb1-analysis-study
isxconcat-sess -sf sessid -analysis sm-gamma-fwhm5 -c odd-v-0 -o group-analysis -hemi lh
}}}

This creates the same files as the volume-based command, except that
they have "lh" or "rh" preprepended instead of "tal".
{{{
ls group-analysis/sm-gamma-fwhm5
ls group-analysis/sm-gamma-fwhm5/odd-v-0
}}}

You should see the following files: lh.meanfunc.nii, lh.mask.nii,
lh.fsnr.nii, lh.ces.nii, lh.cesvar.nii.

=== Group GLM Analysis (Surface) ===

The commands for the GLM analysis are the same as for the volume-based
except that "--surf fsaverage lh" is added to the command line:

Cd into the directory with the assembled data:

{{{
cd $FSFTUTDIR/fb1-analysis-study/group-analysis-tut/sm-gamma-fwhm5/odd-v-0
}}}

==== Random Effects (RFx, OLS) (Surface) ====
{{{
mri_glmfit --y lh.ces.nii --osgm --mask ../lh.mask.nii --glmdir lh.rfx.osgm --nii \
  --surf fsaverage lh
}}}

==== Weighted Random Effects (WRFx, WLS) (Surface) ====
{{{
mri_glmfit --y lh.ces.nii --osgm --glmdir lh.wrfx.osgm --nii --mask ../lh.mask.nii \
    --wls lh.cesvar.nii --surf fsaverage lh
}}}

==== Fixed Effects (FFx) (Surface) ====
{{{
mri_glmfit --y lh.ces.nii --osgm --glmdir lh.ffx.osgm --nii --mask ../lh.mask.nii \
   --yffxvar lh.cesvar.nii --ffxdofdat ../ffxdof.dat --surf fsaverage lh
}}}

=== Output and visualization (Surface) ===

For convenience, we concatenate the different analyses into a single
file, then visualize with tksurfer:

{{{
mri_concat lh.rfx.osgm/osgm/sig.nii lh.wrfx.osgm/osgm/sig.nii lh.ffx.osgm/osgm/sig.nii \
        --o lh.all.sig.nii
tksurfer fsaverage lh inflated -annot aparc -overlay lh.all.sig.nii -fthresh 2
}}}


=== Correction for Multiple Comparisons/Cluster Analysis (Surface) ===

In principle, correction for multiple comparisons on the surface is
done in the same way as volume-based, but in implementation is
considerably more difficult. We have several methods for correcting on
the surface, including False Discovery Rate (FDR), Monte-Carlo
Simulation, and Permutation Simulation. The Guassian Random Field
(GRF) is still being developed (as are the tutorials for it).


= Play with Your Data =

In this section, you are encouraged to experiment with the data. These
analyses have not already been pre-computed, but they should not take
long to run.

== Change the smoothing level ==

Use this series of commands to analyze the data smoothed by 10 mm
FWHM. The preproc-sess command will create fmcsm10.nii. In
mkanalysis-sess, use the -clone option to copy the parameters from an
existing analysis, then just change the Func Stem to fmcsm10.

{{{
cd $FSFUTDIR/fb1-analysis-study
preproc-sess -s mgh-101.1 -fwhm 10
mkanalysis-sess -gui -clone sm-gamma-fwhm5 -analysis sm-gamma-fwhm10
selxavg3-sess -s mgh-101.1 -a sm-gamma-fwhm10
tkmedit-sess -s mgh-101.1 -analysis sm-gamma-fwhm10 \
  -c odd-v-0 -c even-v-0 -c odd+even -c odd-v-even -aparc+aseg

}}}

== Change noise/nuisance model ==

Use this series of commands to analyze the data with little
compensation for noise and nuisance variables. Change the analysis by
turning off whitening, set the Polynomial Order to 0, and uncheck the
MC Regressors.

{{{
cd $FSFUTDIR/fb1-analysis-study
mkanalysis-sess -gui -clone sm-gamma-fwhm5 -analysis sm-gamma-fwhm5B
selxavg3-sess -s mgh-101.1 -a sm-gamma-fwhm5B
tkmedit-sess -s mgh-101.1 -analysis sm-gamma-fwhm5B \
  -c odd-v-0 -c even-v-0 -c odd+even -c odd-v-even -aparc+aseg

}}}
This gives the vector motion at each time point for each run. Note
that it is always positive because this is a magnitude. It is also 0
at the middle time point because the middle time point is used as the
reference.

There are no rules for how much motion is too much motion. Generally
speaking, suddent motions are the worst as are task-related motion.

=== Functional-Anatomical Cross-modal Registration ===

You can get a summary of registration quality using the following command:
{{{
tkregister-sess -s sess01 -fsd bold -per-run -bbr-sum
}}}

This prints out a number of each run found in the
register.dof6.mincost mentioned above. This will be a number between 0
and 1, with 0 being perfect and 1 being terrible. Actual values depend
upon exactly how you have acquired your data. Generally, anything over
0.8 indicates that something is probably wrong.

You can view problematic registrations using the following command:
{{{
tkregister-sess -s sess01 -fsd bold -per-run
}}}

== Notes on Preprocessing ==

* Any command that takes a '-s sessid' argument will take a '-sf sessidfile' argument


= Analysis =

top | previous

FreeSurfer Slides

  1. Functional Analysis with FS-FAST

This tutorial steps you through the analysis of an fMRI data set with the FreeSurfer Functional Analysis Stream (FSFAST) version 5.1, from organizing the data to group analysis.



Tutorial Data Description

The data being analyzed were collected as part of the Functional Biomedical Research Network (fBIRN, www.nbirn.net).

  • Working-memory paradigm with distractors
  • 18 subjects
  • Each subject has 1 run (except sess01 which has 4 runs)
  • Collected at MGH Bay 4 (3T Siemens)
  • FreeSurfer anatomical analyses



Getting and Organizing the Tutorial

If you do not have the tutorial data set up, then consult the FsFastTutorialData page. You will need to set the FSFTUTDIR environment variable. NOTE: if you are taking a class at MGH, the data have already been set up on your computer.

cd into the tutorial data directory and run ls:

cd $FSFTUTDIR

This is the Project Directory. You will run most of the FSFAST commands from the Project Directory.

Run ls to see the contents of the Project Directory.

ls

You will see 18 folders with names like "sess09". These are the 18 subjects. There are some other files and folders there but don't worry about them right now.

Understanding the FS-FAST Directory Structure

All of these sessions have been analyzed with the exception of sess01.noproc. This session is what the directory structure should look like immediately prior to beginning analysis. This includes:

  • Directory stucture
  • Raw data
  • subjectname file
  • Parasigm files for each run

The directory structure and raw data are usually created by "unpacking" the data with dcmunpack or unpacksdcmdir, but it could also be done by hand. The subjectname file and paradigm files must be added manually.

The 'Session' Folder

The folder/directory where all the data for a session are stored is called the 'session' or the 'sessid'. There may be more than one session for a given subject (eg, in a longitudinal analysis). Go into the sess01.noproc folder and run 'ls':

cd sess01.noproc
ls

You will see see two folders ('bold' and 'rest') and a file called 'subjectname'.

There is a file called 'sessidlist' that has a list of all the sessions in this project. These types of files are called 'sessid files'. They are convenient for

* Processing multiple sessions with one command-line * Grouping sessions together * Running analyses in parallel

The 'subjectname' File

subjectname is a text file with the name of the FreeSurfer subject as found in $SUBJECTS_DIR (ie, the location of the anatomical analysis).

View the contents of the text file (with 'cat', 'more', 'less', 'gedit', or 'emacs') and verify that this subject is in the $SUBJECTS_DIR.

NOTE: it is important that the anatomical data and the functional data be from the same subject. The contents of the subjectname file is the only link! Make sure that it is right! There is a check for this below.

Functional Subdirectories (FSDs)

The other two directories (bold and rest) are 'functional subdirectories' (FSDs) and contain functional data. If you 'ls rest' you will see '001'. If you 'ls bold' you will see '001 002 003 004'. Each of these is a 'run' of fMRI data, ie, all the data collected from a start and stop of the scanner.

Raw Data

Go into the first run of the bold directory with 'cd bold/001' and run ls. You will see 'f.nii.gz wmfir.par workmem.par'. The raw data is stored in f.nii.gz (compressed NIFTI) and is directly converted from the DICOM file. Examine this file with mri_info:

mri_info --dim f.nii.gz
mri_info --res f.nii.gz

The first command results in '64 64 30 142'. This is the dimension of the functional data. Since it is functional, it has 4 dimensions: 3 spatial and one temporal (ie, 64 rows, 64 cols, 30 slices, and 142 time points or TRs or frames). The second command results in '3.438 3.437 5.000 2000.000'. This is the resolution of the data, ie, each voxel is 3.438mm by 3.437mm by 5.000mm and the TR is 2000ms (2 sec).

View the functional data with:

tkmedit -f f.nii.gz -t f.nii.gz

Click on a point to view the waveform at that point.

Paradigm Files

The workmem.par and wmfir.par files are paradigm files. They are text files that you create that indicate the stimulus schedule (ie, which stimulus was presented when).

View workmem.par. Each row indicates a stimulus presentation. You will see that each row has 5 columns. The columns are:

  • Stimulus Onset Time (sec)
  • Numeric Stimulus Identifier
  • Stimulus Duration (sec)
  • Weight (usually 1)
  • Text Stimulus Identifier (redundant with Numeric Stimulus Identifier)

The Stimulus Onset Time is the onset relative to the acquiistion time of first time point in f.nii.gz. The Numeric and Text Stimulus Identifiers indicate which stimulus was presented. The Stimulus Duration is the amount of time the stimulus was presented.

In this case, there are 5 event types:

  • Encode - encode phase
  • EDistrator - emotional distractor
  • EProbe - probe after emotional distractor
  • NDistrator - neutral distractor
  • NProbe - probe after neutral distractor

Note two things: (1) Not all the time is taken up, and (2) Basline/Fixation is not explicitly represented. By default, any time not covered by stimuluation is assumed to be baseline.

  • What time was the third Encode presented?
  • What time was the last Neutral Distractor presented?
  • Is this timing for all runs the same?

Preprocessing

Once the data have been arranged in the above directory structure and naming convention, they are ready to be preprocessed. In FS-FAST, it is assumed that each data set will be analyzed in three ways:

  • Left Cortical Hemisphere
  • Right Cortical Hemisphere
  • Subcortical Structures

You will need to decide how much to smooth the data and whether you want to do slice-timing correction. In this analysis, we will smooth the data by 5mm Full-Width/Half-Max (FWHM) and correct for slice timing. The slice-timing for this particular data set was 'Ascending', meaning that the first slice was acquired first, the second slice was acquired second, etc. To preprocess the data, run:

preproc-sess -s sess01 -surface fsaverage lhrh -mni305 \
   -fwhm 5 -stc up -fsd bold -per-run 

This data has already been preprocessed, so it should just verify that it is up-to-date and return. This command has several arguments:

  • -s sess01 : indicates which session to preprocess
  • -surface fsaverage lhrh : indicates that data should be sampled to the left and right hemispheres of the 'fsaverage' subject
  • -mni305 : indicates that data should be sampled to the mni305 volume (2mm isotropic voxel size)
  • -fwhm 5 : smooth by 5mm FWHM
  • -stc up : slice-timing correction with ascending ('up') slice order
  • -fsd bold : indicate the functional subdirectory
  • -per-run : register each run separately

This command does a lot (and it can take quite a long time to run). To understand what it does, we will go back into one of the run directories and see what it creates. To do this, 'cd sess01/bold/001' and type 'ls'. This directory previously held only f.nii.gz, workmem.par, and wmfir.par; now there are a lot of files, each indicative of a different preprocessing stage. Now type 'ls -ltr'. This command sorts the files by creation time with the oldest at the top and the newest at the bottom. The preprocessing is progressive, meaning that the output of one stage is the input to the next.

Template

This stage creates template.nii.gz (and template.log). This is the middle time point from the raw functional data (f.nii.gz). This is the reference used to motion correct and register the functionals for this run to the anatomical. It is also used to create masks of the brain.

* Verify that it has the same dimension and resolution as the raw data using mri_info.

Masking

The masks for this run are stored in the 'masks' directory. Run 'ls -ltr masks'. You will see a file called 'brain.nii.gz'. This is a binary mask created using the FSL BET program. There is also a file called 'brain.e3.nii.gz' which is the mask eroded by three voxels. These have the same dimensions as the template. View the masks with:

tkmedit -f template.nii.gz -overlay masks/brain.nii.gz -fthresh 0.5
tkmedit -f template.nii.gz -overlay masks/brain.e3.nii.gz -fthresh 0.5

The brain.nii.gz is used to constrain voxel-wise operations. The eroded mask (brain.e3.nii.gz) is used to compute the mean functional value used for intensity normalization and global mean time course. There are other masks there that we will get to later.

Intensity Normalization and Global Mean Time Course

By default, FSFAST will scale the intensities of all voxels and time points to help assure that they are of the same value across runs, sessions, and subjects. It does this by dividing by the mean across all voxels and time points inside the brain.e3.nii.gz mask, then multiplying by 100. This value is stored in global.meanval.dat. This is a simple text file which you can view. At this point, this value is stored and used later. A waveform is also constructed of the mean at each time point (text file global.waveform.dat). This can be used as a nuisance regressor.

* What are the global means for runs 1, 2, 3, and 4?

Functional-Anatomical Cross-modal Registration

The next six files (init.register.dof6.dat, register.dof6.dat, register.dof6.dat.mincost, register.dof6.dat.sum, register.dof6.dat.log, register.dof6.dat.param) deal with the registration from the functional to the same-subject FreeSurfer anatomical. There are only two files here that are really important: register.dof6.dat and register.dof6.dat.mincost.

  • register.dof6.dat - is a text file that contains the registration matrix.
  • register.dof6.mincost - is a text file that contains a measure of
    • the quality of the registration.

The registration is will be revisited below when we talk about Quality Assurance

Motion Correction (MC)

The motion correction stage produces these files: fmcpr.mat.aff12.1D, fmcpr.nii.gz, mcprextreg, mcdat2extreg.log, fmcpr.nii.gz.mclog, fmcpr.mcdat. There are only three important file here:

  • fmcpr.nii.gz -- this is the motion corrected functional data. It has the same size and dimension as the raw functional data.
  • fmcpr.mcdat - text file of the amount of motion at each time point. This is important for Quality Assurance (see below).
  • mcprextreg - text file of the motion correction parameters assembled into an orthogonalized matrix that can be used as nuisance regressors
  • Verify that fmcpr.nii.gz has the same dimension as f.nii.gz using mri_info.

Slice-Timing Correction (STC)

Slice-timing corretion compensates for the fact that each of the 30 slices was acquired separately over the course of 2 sec. It does this by interpolating between time points to align each slice to the time of the middle of the TR. The file created with this is fmcpr.up.nii.gz (and fmcpr.up.nii.gz.log).

  • Verify that fmcpr.up.nii.gz has the same dimension as f.nii.gz using mri_info.

Resampling to Common Spaces and Spatial Smoothing

At this point, the functional data has stayed in the 'native functional space', ie, 64x64x30, 3.4x3.4x5mm3. Now it will be sampled into the 'Common Space'. The Common Space is a geometry where all subjects are in voxel-for-voxel registration. There are three such spaces in FSFAST:

  • Left hemisphere of fsaverage (fmcpr.up.sm5.fsaverage.lh.nii.gz)
  • Right hemisphere of fsaverage (fmcpr.up.sm5.fsaverage.rh.nii.gz)
  • Volume of fsaverage (MNI305 space) - for subcortical analyses (fmcpr.sm5.mni305.2mm.nii.gz)

Each of these is the entire 4D functional data set resampled into the common space. The spatial smoothing is performed after resampling. Surface-based (2D) smoothing is used for the surfaces; 3D for the volumes.

Check the dimensions of the MNI305 space volume:

mri_info --dim fmcpr.up.sm5.mni305.2mm.nii.gz
mri_info --res fmcpr.up.sm5.mni305.2mm.nii.gz

The dimension will be '76 76 93 142' meaning that there are 76 columns, 76 rows, 93 slices but still 142 time points (same as the raw data). The resolution will be '2.0 2.0 2.0 2000' meaning that each voxel is 2mm in size and the TR is still 2 sec. The transformation to this space is based on the 12DOF talairach.xfm created during the FreeSurfer reconstruction.

Check the dimensions of the Left Hemisphere 'volume':

mri_info --dim fmcpr.up.sm5.fsaverage.lh.nii.gz
mri_info --res fmcpr.up.sm5.fsaverage.lh.nii.gz

The dimension is '163842 1 1 142'. This 'volume' has 163842 'columns', 1 'row', and 1 'slice' (still 142 time points). You are probably confused right now. That's ok, it's natural. At this point the notion of a 'volume' has been lost. Each 'voxel' is actually a vertex (of which there are 163842 in the left hemisphere of fsaverage). Storing it in a NIFTI 'volume' is just a convenience.

The 'resolution' is '1.0 1.0 1.0 2000'. The values for the first 3 dimensions are meaningless because there are no columns, rows, or slices on the surface so the distances between them are meaningless. The last value indicates the time between frames and is still accurate (2 sec).

The transformation to this space is based on the surface-based intersubject registration created during the FreeSurfer reconstruction.

Quality Assurance

Motion Correction

The motion plots can be viewed with:

plot-twf-sess -s sess01 -fsd bold -mc 

This gives the vector motion at each time point for each run. Note that it is always positive because this is a magnitude. It is also 0 at the middle time point because the middle time point is used as the reference.

There are no rules for how much motion is too much motion. Generally speaking, suddent motions are the worst as are task-related motion.

Functional-Anatomical Cross-modal Registration

You can get a summary of registration quality using the following command:

tkregister-sess -s sess01 -fsd bold -per-run -bbr-sum

This prints out a number of each run found in the register.dof6.mincost mentioned above. This will be a number between 0 and 1, with 0 being perfect and 1 being terrible. Actual values depend upon exactly how you have acquired your data. Generally, anything over 0.8 indicates that something is probably wrong.

You can view problematic registrations using the following command:

tkregister-sess -s sess01 -fsd bold -per-run

Notes on Preprocessing

* Any command that takes a '-s sessid' argument will take a '-sf sessidfile' argument

Analysis

FsFastTutorialV5.1 (last edited 2016-09-01 16:35:55 by AllisonMoreau)