Detrending FreeSurfer Data: Difference between revisions

From CCN Wiki
Jump to navigation Jump to search
Line 3: Line 3:
A script has been written called detrend.sh that removes the linear trend in your BOLD data. The latest version of this script can be found in /usr/local/sbin (which should be in your $PATH):
A script has been written called detrend.sh that removes the linear trend in your BOLD data. The latest version of this script can be found in /usr/local/sbin (which should be in your $PATH):
==== detrend.sh ====
==== detrend.sh ====
#!/bin/bash
  #!/bin/bash
  USAGE="Usage: detrend.sh filepattern surface sub1 ... subN"
  USAGE="Usage: detrend.sh filepattern surf sub1 ... subN"
  if [ "$#" == "0" ]; then
  if [ "$#" == "0" ]; then
echo "$USAGE"
        echo "$USAGE"
exit 1
        exit 1
  fi
  fi
 
  #first parameter is the filepattern for the .nii.gz time series to be detrended, up to the surface indicator
  #first parameter is the filepattern for the .nii.gz time series to be detrended, up to the surface indicator
  #e.g., fmcpr.sm6.$surf.?h.nii.gz would use fmcpr.sm6 as the filepattern
  #e.g., fmcpr.siemens.sm6.fsaverage.?h.nii.gz would use fmcpr.siemens.sm6 as the filepattern, fsaverage as the surf
  filepat="$1"
  filepat="$1"
  shift
  shift
   
   
  #second parameter should be specified either as '''self''' or '''fsaverage'''
  #second parameter should be specified either as self or fsaverage
  surf="$1"
  surf="$1"
  shift
  shift
   
   
  #after the shift command, all the arguments are shifted down one place and the first two arguments
  #after the shift command, all the arguments are shifted down one place and the first argument (the filepattern)  
#(the filepattern and surface)  
  #falls off the list. The remaining arguments should be subject_ids
  #fall off the list. The remaining arguments should be subject_ids
  subs=( "$@" );
  subs=( "$@" );
  hemis=( "lh" "rh" );
  hemis=( "lh" "rh" );
Line 32: Line 32:
         echo "${source_dir} does not exist!"
         echo "${source_dir} does not exist!"
     else
     else
cd ${source_dir}
        cd ${source_dir}
readarray -t runs < runs
        readarray -t runs < runs
for r in "${runs[@]}"; do
        for r in "${runs[@]}"; do
if [ -n "${r}" ]; then
                if [ -n "${r}" ]; then
#the -n test makes sure that the run number is not an empty string
                #the -n test makes sure that the run number is not an empty string
#caused by a trailing newline in the runs file
                #caused by a trailing newline in the runs file
for hemi in "${hemis[@]}"; do
                        for hemi in "${hemis[@]}"; do
cd ${source_dir}/${r}
                                cd ${source_dir}/${r}
pwd
                                pwd
        #subject_id does exist. Detrend
                                #subject_id does exist. Detrend
if [ "${surf}" = "self" ]; then
                                if [ "${surf}" = "self" ]; then
SURFTOUSE=${sub}
else
SURFTOUSE=fsaverage
fi
                mri_glmfit --y ${source_dir}/${r}/${filepat}.${surf}.${hemi}.nii.gz \
        --glmdir ${source_dir}/${r}/${hemi}.detrend \
--qa --save-yhat --eres-save \
--surf ${SURFTOUSE} ${hemi}
mv ${source_dir}/${r}/${hemi}.detrend/eres.mgh ${source_dir}/${r}/${filepat}.${hemi}.mgh
        done
   
   
                      <span style="color:red;">#NEW! Now the script also detrends the mni305 file
                                        SURFTOUSE=${sub}
                                else
                                        SURFTOUSE=fsaverage
                                fi
                                mri_glmfit --y ${source_dir}/${r}/${filepat}.${surf}.${hemi}.nii.gz \
                                --glmdir ${source_dir}/${r}/${hemi}.detrend \
                                --qa --save-yhat --eres-save \
                                --surf ${SURFTOUSE} ${hemi}
                                mv ${source_dir}/${r}/${hemi}.detrend/eres.mgh ${source_dir}/${r}/${filepat}.${hemi}.mgh
                        done
                        #now detrend the mni305 file
                         mri_glmfit --y ${source_dir}/${r}/${filepat}.mni305.2mm.nii.gz \
                         mri_glmfit --y ${source_dir}/${r}/${filepat}.mni305.2mm.nii.gz \
                         --glmdir ${source_dir}/${r}/mni305.detrend \
                         --glmdir ${source_dir}/${r}/mni305.detrend \
                         --qa --save-yhat --eres-save
                         --qa --save-yhat --eres-save
                         mv ${source_dir}/${r}/mni305.detrend/eres.mgh ${source_dir}/${r}/${filepat}.mni305.2mm.mgh</span>
                         mv ${source_dir}/${r}/mni305.detrend/eres.mgh ${source_dir}/${r}/${filepat}.mni305.2mm.mgh
fi
                fi
done
        done
     fi
     fi
  done
  done

Revision as of 21:31, 11 June 2019

Over the course of a run, there can be a linear drift in the signal in different regions of the brain. There are many possible causes for this that have nothing to do with any interesting aspect of your data -- in other words, this linear drift is a nuisance artifact. The second step is to remove this signal drift from the data because it can introduce spurious correlations between two unrelated time series. You can see this for yourself in a quick experiment you could whip up in Excel: take two vectors of 100 randomly generated numbers (e.g., randbetween(1,99)). They should be uncorrelated. Now add 1, 2, 3, ... , 99, 100 to the values in each vector. This simulates a linear trend in the data. You shouldn't be surprised to find that the two vectors are now highly and positively correlated!

A script has been written called detrend.sh that removes the linear trend in your BOLD data. The latest version of this script can be found in /usr/local/sbin (which should be in your $PATH):

detrend.sh

 #!/bin/bash
USAGE="Usage: detrend.sh filepattern surf sub1 ... subN"

if [ "$#" == "0" ]; then
       echo "$USAGE"
       exit 1
fi

#first parameter is the filepattern for the .nii.gz time series to be detrended, up to the surface indicator
#e.g., fmcpr.siemens.sm6.fsaverage.?h.nii.gz would use fmcpr.siemens.sm6  as the filepattern, fsaverage as the surf
filepat="$1"
shift

#second parameter should be specified either as self or fsaverage
surf="$1"
shift

#after the shift command, all the arguments are shifted down one place and the first argument (the filepattern) 
#falls off the list. The remaining arguments should be subject_ids
subs=( "$@" );
hemis=( "lh" "rh" );

for sub in "${subs[@]}"; do
   source_dir=${SUBJECTS_DIR}/${sub}/bold
   echo ${source_dir}
   if [ ! -d ${source_dir} ]; then
       #The subject_id does not exist
       echo "${source_dir} does not exist!"
   else
       cd ${source_dir}
       readarray -t runs < runs
       for r in "${runs[@]}"; do
               if [ -n "${r}" ]; then
               #the -n test makes sure that the run number is not an empty string
               #caused by a trailing newline in the runs file
                       for hemi in "${hemis[@]}"; do
                               cd ${source_dir}/${r}
                               pwd
                               #subject_id does exist. Detrend
                               if [ "${surf}" = "self" ]; then

                                       SURFTOUSE=${sub}
                               else
                                       SURFTOUSE=fsaverage
                               fi
                               mri_glmfit --y ${source_dir}/${r}/${filepat}.${surf}.${hemi}.nii.gz \
                               --glmdir ${source_dir}/${r}/${hemi}.detrend \
                               --qa --save-yhat --eres-save \
                               --surf ${SURFTOUSE} ${hemi}
                               mv ${source_dir}/${r}/${hemi}.detrend/eres.mgh ${source_dir}/${r}/${filepat}.${hemi}.mgh
                       done
                       #now detrend the mni305 file
                       mri_glmfit --y ${source_dir}/${r}/${filepat}.mni305.2mm.nii.gz \
                       --glmdir ${source_dir}/${r}/mni305.detrend \
                       --qa --save-yhat --eres-save
                       mv ${source_dir}/${r}/mni305.detrend/eres.mgh ${source_dir}/${r}/${filepat}.mni305.2mm.mgh
               fi
       done
   fi
done

Running the Script

Before running this script, you will need to create a text file called 'runs' in the bold/ directory for each subject's dataset, e.g.,

  • FS_T1_501/
    • bold/
      • runs
      • 005/
      • 006/

The runs file simply lists each run folder on its own line:

005
006

The detrend.sh script uses this file to determine the folders containing the data to be detrended. If this file doesn't already exist, you can manually generate it in any text editor (e.g., nano runs or gedit runs), but the quickest method takes advantage of the fact that the run folders all start with 0, and uses common command-line utilities:

SUBJECT=FS_501
cd ${SUBJECTS_DIR}/${SUBJECT}/bold
#if there are fewer than 10 runs of BOLD data, then the run directories will probably have 2 leading zeros
ls -1 | grep "^00*" > runs

Assuming all your subject folders have the same run folders to detrend, you would detrend multiple subjects using detrend.sh, specifying a file pattern for the source data (i.e., the name of the preprocessed files generated by FS-FAST, omitting anything after the '?h' hemisphere identifier), followed by a list of subject IDs:

#A SPECIFIC EXAMPLE: (note these parameters may differ substantially from what you would be typing in)
detrend.sh fmcpr.sm6 self FS_T1_501 FS_T2_501 FS_T1_505 FS_T2_505
#For a more generalizable example of how you should call this function, see the section below using variables

The gist is that it calls the mri_glmfit function and saves the residuals after the linear trend has been removed from the data. Multiple files are generated in ?h.detrend/ directories in each run directory. The detrended data is subsequently copied back to the run directory as a new file called ${filepat}.?h.mgh, where ${filepat} is whatever file pattern you provided to the script (note that the source data are .nii.gz files, whereas the detrended data are .mgh files).

Using variables

As I mentioned in the comments in the sample code snippet above, what you type depends entirely on the filenames, which in turn depend entirely on how the data were preprocessed. You can use environment variables to help walk you through figuring out the correct file patterns. Also, a handy shortcut exists if you happen to have a subjects file in $SUBJECTS_DIR. Putting these two techniques together:

FILEPATTERN=fmcpr #the preprocessed files will almost always be called fmcpr
SMOOTHING=".sm4" #how much smoothing did you use when you ran preproc-sess?
SLICETIME=".up" #[".up" | ".down" | ".siemens" | OMIT ]
SURFACE=fsaverage #[self | fsaverage]

detrend.sh ${FILEPATTERN}${SLICETIME}${SMOOTHING} $SURFACE `cat subjects`

This will execute the detrend.sh script on all the subjects listed in the subjects text file, using the fsaverage surface. The assemblage of variables will look for files named fmcpr.up.sm4.*