Region of Interest Analysis with FSL
  Introduction

This page describes how to analyze fMRI data from a single individual using FSL. FSL is a powerful and free tool to measure task-related changes in the brain's blood flow. This tutorial also shows how to conduct region of interest analysis using Matlab code written by Russell Poldrack. To complete this tutorial you will need:

  • To conduct FSL single subject data processing you will need
    1. FSL installed on your computer (you must first install Cygwin if you are running Windows).
    2. The sample dataset. This includes a fMRI session, and text file that describes the events that occurred during the fMRI session.
  • To complete the Peristimulus Plot for a Region of Intererest you will need all of the above as well as:
    1. Matlab must be installed.
    2. Russell Poldrack's Matlab scripts must be in Matlab's path. We have made a number of modifications to allow these scripts to run under Windows as well as Linux.
    3. The sed stream editor must be installed. To see if you have this, type 'which sed' from the Unix command line (Unix) or the Cygwin command line (Windows). If sed is available, you should see a response like this '/usr/bin/sed'.

About the sample dataset

This was a simple finger tapping task. For periods of 18 seconds, the screen rapidly flashed the word 'rest' on the screen, and the participant was asked to do nothing, next the participant saw the words 'R1' and 'R2' flash on the screen rapidly for 18 seconds, and was expected to press the right thumb and index finger. These 36-second rest-activation cycles repeated for the entire scanning session (about 4 minutes). The first 5 fMRI scans were discarded (as these appear brighter as the T1 effects are not saturated), with 130 3D volumes saved to disk (one volume every 1.78 seconds). The first finger task began 11.3 seconds after the first saved fMRI dataset. The image below shows the task (blue square wave) and a measured change in image intensity. The lower panel shows how an average response can be generated, showing the mean change in signal in response to a stimulus. In particular, in this tutorial we will generate a peristimulus plot using a pre-defined region of interest that marks the location of the motor hand area.


Above: the motor hand area.

FSL single subject data processing

Steps for functional MRI analysis using FSL version 3.1 running under Windows XP:August, 19th, 2004.

  1. Start FSL by clicking on the fsl shortcut within the computer desktop, or by starting typing 'fsl' from the command prompt (or from the CygWin command prompt if you are running Wndows). Note that if you type 'fsl &', the command console under which FSL is running will not be "frozen" by FSL, but can be used for other commands.
  2. By typing fsl, a user interface (the one with a fossil fish on its top) appears. This interface contains several options, but the one to be used now for fmri is 'FEAT FMRI Analysis' button. Press the corresponding button and this will open the Feat interface.(The forthcoming steps are described for use in FEAT version 5.1)
  3. The feat interface has two top buttons. The following analysis is a first level analysis (within one single subject), and it is a full analysis (involving preprocessing of the data, the statistical spinnings of fmri analsysi and post-stats registrations).
  4. The FEAT interface five tabs across the top portion of the screen, oriented in a horizontal line. The fmri analysis can be set by defining the parameters within each one of those buttons, in a left to right order:
    • within Misc Tab: - leave 'Featwatcher' checked to gauge the progress of the fmri analysis. - set the delay to start the analysis to a number different to zero: we want to compute this analysis immediately.
    • within Data Tab: - Number of analysis is set to one. - Press select 4D data. Browse to the folder in which the 4D data file is stored and select the .hdr of the corresponding file *See footnote. After selecting the 4D file, the number of volumes will automatically be displayed. In our case 130. - enter the TR. In our case 1.77939. - it is not necessary to delete any volumes. - set the high pass filter cutoff to 128 s. - You can set the output directory. If you leave it blank, it will write the results in the same folder where the 4D data is. FEAT write the results in a folder called "4Ddataname.feat". In our case fHAM26S.feat. If you re-run the analysis, it will create a second folder called fHAM26S+.feat, and so on.
    • within Pre-stats: - Slice timing correction should be set to regular up (0,1,2...n-1). - Motin correction is done with Motion Correction using FMRIB's Linear Image Registration Tool - MCFLIRT. - Maintain Brain Extraction checked. - Spatial smoothing FWHM set to 6mm. - Intensity Normalization (equivalent to Grand Mean Scaling) unchecked. - Temporal filtering Highpass checked.
    • Within Stats Tab - Simple model setup could be used for alternating designs like ABAB, or ABCABC. Even though our design is an alternating block design, to ensure accuracy, press the 'Full Model Setup' button (that can also be used for event related designs). - Clicking full model setup opens a General Linear model interface. Within this interface, follow the steps: - The number of original EV (original EVS - explanatory variables, or simply put: conditions) is set to 1, because we are interested in the motor actication of the "action" condition. - Each EV is setup separately. - The basic shape of the wave form that describes the stimulus that we wish to model is defined by a custom file expliciting the time in seconds of the onset of the stimulus, how much time (also in seconds) the stimulus last, and the value of the input. This requires a .txt file containing a matrix under which the columns are the variables stated above (i.e., first one time of onset, second one the duration, and the thrid one the value) and the rows represent the stimuli
      Contents of cond1_onset_duration.txt
      11.32366 18 1
      47.32366 18 1
      83.32366 18 1
      119.32366 18 1
      155.32366 18 1
      191.32366 18 1
      227.32366 18 1

      This is set as follows: first choose the 'Custom (3 column format)' from the 'Basic Shape' pull-down menu, second browse to find the file cond1_onset_duration.txt *See Footnote, third set the convolution (which is the form of the Hemodynamic response function that will be applied to the basic waveform) to 'Gamma', fourth maintain Phase, Stddev and Mean Lag as defaults (respectively 0,3 and 6 secs), fifth apply temporal filtering and add temporal derivative. Once these parameters have been entered click on 'Contrasts & F-tests' tab. Since we are intersted in the effects of activation of only one condition, contrasts are to be set to one, F-tests to zero, leave OC checked, add the name activation to OC1 and maintain EV1 as one Click view design andn you should see the 2D image of the design matrix. It contains two columns showing different timecourses, with time going down vertically . The first one represents the time course associated with our EV1, and the second one the model to be fit. You can now press 'Done' to close the 'General Linear Model' window.

    • Within Post_stats Tab- Do not use a pre-threshold masking - Maintain defaults for thresholding (Cluster analysis, Z threshold of 2,3, cluster P of 0.01) - Maintain defaults for rendering and contrast masking.
    • Within Registration Tab: check the 'Main structural image' checkbox and insert the correct filename (e.g. C:/sensefsl/aHAM.hdr'). Note: FLIRT works much better if you scalp strip your anatomical scans prior to registration. The sample image has been scalp-stripped using the outstanding BET tool included with FSL. *See Footnote,.
  5. Press the 'Go' button on the bottom left. and Featwatcher will be opened and will show you the progress of the analysis. After the analysis finishes, Feat watcher will bring you a message stating the time it finished and where the folder with the results is.
  6. Look for the folder called fHAM26S.feat that will probably be in the same folder where the 4D data is.
  7. To assess the analysis click on report.html and the internet browser used as the default in your computer will bring the page containing
    • Click the link to view the motion correction report
    • Click the Thresholded activation images: you are redirected to a link containing the cluster list with the corresponding local maxima eiher in voxel space version or in Talairach space version.
    • The time series plot for the whole experiment
    • A summary of the analysis methods and its referencese) the design and covariance matrices

* Footnote: FSL 3.1 working under Windows has an annoying bug, which is to duplicate the name of the path. If you want it to be C:/sensefsl/fHAM26S.hdrit will automatically write C:/sensefsl/C:/sensefsl/fHAM26S.hdr To overcome this problem, just delete the unwanted part.

Peristimulus Plot for a Region of Intererest:

To conduct this next stage you need to have completed the previous steps. You also need to have Matlab installed and must make sure that the relevant .m files are in your Matlab path. We want to plot the activation in the right motor hand region. In particular, we will apply a region of interest that includes the central locations reported for the motor hand area. To run this script, we simply need to edit the Matlab script 'batch_psplot'. This script simply allows you to specify the location of your files and runs two of Russell Poldrack's functions for creating peristimulus plots.

  1. selective_avg : this script will compute average timecourses for each voxel and for each condition. It requires a few settings (parameters). Selective_avg creates a series of images: for each condition a 4D image is created, with each 3D volume in the series representing a different point in time relative to the condition onset.
    1. params.featdir: the location of the FSL analysis you conducted in the previous portion of this tutorial
    2. params.sed_converter: the location of a file that gives instructions for converting a FSL format description of the study into a Matlab-readable version. Note that your Unix/Cygwin system will have to have the program 'sed' installed: for example SUSE Linux does not automatically install SED, so you will need to add it using YaST2.
    3. params.clusterInStandardSpace: [optional]: Selective_avg creates a timeseries output in the original space of the individual's brain (in the .feat\roi) folder. However, if you specifiy a image in standard space, then selective_avg will also create a spatially normalized timeseries (in the .feat\reg\roi folder). This is useful if you want to apply a region of interest that is in standard stereotaxic space. This allows you to use predefined regions of interest as well as allowing you to apply a region of interest for each and every individual in a group analysis. Important: with FSL 3.1 this region of interest must have exactly the same spatial dimensions used for the 'registration' step of image processing, e.g. looking at the previous section you will notice that we coregistered our images to the 'avg152T1_brain' image: an image that is 91x109x91 voxels in size, with between-voxel spacing of 2mm in each dimension. To turn off this function, simply set this value to '';
    4. params.NegWindow: The amount of time prior to stimulus onset that you want to plot. For example, a value of 4 means that you want the software to plot events occurring 4 seconds prior to stimulus onset.
    5. params.PosWindow: Amount of time after stimulus onset that you want to plot. Typically, this will be only particularly meaningful for the first 8-16 seconds for event related designs (as later time-bins will show increasing interference from other stimuli). However, for block designs you probably want to plot beyond the end of the block length to see the full shape of the hemodynamic response.
    6. params.ter: The size (in seconds) of each timebin in the output. This is typically set to the TR (time to collect a 3D volume of data).
  2. extract_clust_ts: This simple script takes the4D timecourse output from selective_avg and maps the average activity in a 3D region of interest. The 3D images that compose the 3D image and the 3D region of interest must have the same spatial orientation and dimensions. The function returns the values ts, which is simply the timeseries for region of interest.
    1. roi_img: a 3D image that highlights a cluster of voxels that compose your region of interest.
    2. cond_avg_img: 4D output from selective_avg.

On my computer, the script looks like this, though you may have to adjust the file locations (e.g. on Unix systems, the data might be in /home/chris/sensefsl/fHAM26S.feat instead of on the C:\ drive). Once you have altered this script and saved the edited file, simply launch Matlab and run 'batch_psplot' (Unix users should use the file 'batch_psplot_linux': this uses a different .sed file that uses Unix line endings).

Contents of batch_psplot.m
%batch_psplot : compute peri-stimulus timecourse for FEAT fMRI data
%change values in batch_psplot.m to match your setup
%first: assign parameters used by selective_avg
params.featdir = 'C:\sensefsl\fHAM26S.feat';
params.sed_converter = 'C:\MATLAB6p5\fsl\designconvert.sed';
params.ClusterInStandardSpace = 'C:\sensefsl\right_handarea.img';
% - note: optional: set to '' to turn off. Provides reference image
% for applyxfm4D (so coregistered output will have correct bounding box)
params.NegWindow = 4; %from 4 seconds prior to onset
params.PosWindow = 30;
%until 30 seconds after stimulus onset
params.ter = 1.77939;
%run selective_avg with parameters: output: one 4D image per condition
%each slice of 4D image repesents different averaged time-point
%output saved in folder 'ROI' within the .FEAT' directory
selective_avg(params);
%optional: conduct cluster (region-of-interest) analysis
%requires
% cond_avg_img: output for a condition of sective_avg
% roi_img: cluster image with same dimensions as cond_avg_img
% notes:
% 1.) you might want to conduct analysis of multiple conditions or clusters
% 2.) you might want to apply to normalized images or raw MRI scans
cond_avg_img = [params.featdir filesep 'reg' filesep 'roi' filesep 'cond1_onset_duration.img'];
roi_img = params.ClusterInStandardSpace;
ts=extract_clust_ts(cond_avg_img, roi_img)
plot(ts);

The function batch_psplot will take a long time to run: in particular reslicing the timeseries into standard space requires a lot of time. Matlab should report the mean intensity for the timeseries, and you can save/copy the values of the timeseries (ts) for group analysis. The values of ts are also graphically shown as a plot (see picture below). Note that the image does not look exactly like the predicted HRF response shown in the figure. The reason for this may be inhibition: at the end of the block there is strong activity in the motor cortex to stop the regular responses.

23 Aug 2004: Leonardo Bonilha and Chris Rorden