ASL/BOLD fMRI Processing

ASL/BOLD fMRI data processing was carried out using the following series of modules:

  • Motion Correction
  • fMRI Preprocessing (Spatial and Temporal filtering)
  • Outlier Identification
  • First-level general linear model
  • Registration of functional volumes to anatomical space
  • Transformation of first-level results to standard space

A detailed description of each of these stages is as follows:

Details of the ASL/BOLD data acquistion [Schmithorst2014] parameters are here: protocols

various Feat .fsf files, etc relating to the ASL/BOLD Stories and Sentences tasks can be found within the cmind-py repository in /cmind/study_specific/ASLBOLD/

Motion Correction

Motion correction is performed using MCFLIRT [Jenkinson2002]. Rather than using a fixed reference frame, multiple calls to mcflirt are made as follows:

  • 1.) Fast initial registration to an initial frame without sinc interpolation to find the frame closest to the mean position
  • 2.) Rerun mcflirt using this “optimal” reference frame and with sinc interpolation to reduce blurring
  • 3.) Take a timeseries average of the reference frames and then run mcflirt a 3rd time using this high SNR timeseries average as the reference

For ASL/BOLD data, there are two sets of 4D timeseries data. In this case, stage 1 is performed independently for the ASL and BOLD timeseries. The average of the motion estimates from both the ASL and BOLD data is used to estimate the frame closest to the mean position. This common reference frame is then used in stage 2.

The estimated motion parameters from stage 3 are stored for future use as nuisance regressors during first level analysis.

For the ASL data, a mean (tag-control) volume is also formed to allow for quick verification of ASL labeling.

Two sets of global motion parameters were generated by MCFLIRT and stored in the processed tar.gz files and can be queried in the database using the following variable names: absolute motion correction and relative motion correction. Timecourse mean absolute motion relative to the reference frame is stored in the file “<ASL|BOLD>_<Stories|Sentences>_mcf_abs_mean.rms”. The mean relative motion between adjacent timepoints is stored in “<ASL|BOLD>_<Stories|Sentences>_mcf_rel_mean.rms”.

Pipeline Function Details:

fMRI Preprocessing

The preprocessing pipeline is closely matched to the standard fMRI preprocessing stream used in FSL [Smith2004], [Jenkinson2012]. This includes brain extraction via BET [Smith2002] (threshold of 0.3), spatial smoothing to 8 mm via SUSAN [Smith1997], and high-pass temporal filtering with a cutoff of 104 seconds. Each 4D data set is normalised by a single scaling factor (grand mean scaling). The same preprocessing stages are repeated for both the ASL and BOLD timeseries.

Slice timing correction was not performed. This is because the paradigms consist of long blocks (64 second ON, 64 second OFF) and the slices are all acquired within a very short time interval (555 ms for ASL, 1.07 seconds for BOLD) within the TR. Under these conditions, slice timing correction would have very little effect and small shifts in the BOLD signal are already taken into account via the temporal derivative of the task regressor.

The parameters listed above describe those used for the processed data stored in the database. In general it is also possible to pass a feat_params.csv file specifying a different spatial smoothing, etc.

Pipeline Function Details:

Outlier Detection

All first-level data stored in the database used the following two sets of nuisance regressors:

    1. six motion regressors as determined during motion correction (3 rotations, 3 translations)
    1. outlier rejection regressors (corresponding to columns in the design matrix that are zero, with only a single one at the timepoint to be rejected)

Outlier timepoints are identified using the fsl_motion_outliers script from FSL, with slight modification to use the timeseries mean rather than a particular timepoint as the reference volume. The metric used is the root-mean-square (RMS) intensity difference between each 3D volume and a reference volume. This is closely related to the DVARS metric proposed by Power et. al. [Power2012]. Slow trends in RMS intensity difference due to drift are removed by first order differencing of the RMS error metric. Outliers are then identified using 1.5 times the inter-quartile range:

  • threshold = p75 + 1.5 * (p75 - p25)

where p25 and p75 are the 25th and 75th percentiles, respectively of the RMS intensity differences.

Pipeline Function Details:

First-level GLM

General linear model (GLM) based first-level analysis is performed using FILM with local autocorrelatoin correction [Woolrich2001].

ASL data is NOT pairwise subtracted prior to processing, but instead regressors corresponding to tag/control baseline and tag/control modulated versions of the paradigm are generated (see: [Mumford2006]).

The design matrix and contrasts for the ASL timeseries and the Stories task is illustrated below. The regressor in the first column is the ASL baseline. The second column is the BOLD task regressor. The third column is the first derivative of the second (to allow small temporal shifts). The fourth column is the product of mean-centered ASL baseline and BOLD task regressor. This represents the ASL task regressor.


ASL Stories design matrix

The following figure demonstrates the covariance among these regressors, demonstrating that they are nearly independent.


covariance matrix

Note #1: The Sentences task design matrix is identical except the order of task (ON) and rest (OFF) in the task regressors is swapped.

Note #2: For the Siemens data from UCLA, the ASL control/tag order is opposite of that on the Philips sequence implemented at CCHMC. The ASL-related regressors are modified to reflect this when analyzing the UCLA data.

Post-stats thresholding was performed on a voxelwise basis at a familywise error (FEW) rate of 0.05.

Registration to Structural Space

The individual volumes of the ASL/BOLD timeseries have very little contrast, making registration more difficult. For this reason the mean ASL subtraction volume, which has high gray/white contrast was used for registration to the T1-weighted structural volume. For subjects without fieldmap data, this is an affine boundary-based registration [Greve2009] as implemented in FLIRT [Jenkinson2002]. For those subjects with fieldmap data, epi unwarping is also applied via FUGUE [Jenkinson2001b].

Normalization to Standard Space

The warp from subject to standard space is combined with the warp from the structural image to MNI space to give parameter estimate and variance maps (FSL’s copes & varcopes) in standard space. Although ANTs is used for some of these stages, the outputs follow FSL naming conventions so that each subjects first-level results can be used inputs to a standard FSL-based second-level analysis.