Deletions are marked like this. | Additions are marked like this. |
Line 1: | Line 1: |
= This Page gathers information about MRI processing (functional, DTI) performed in the Language and Reading Research Lab (LRRL) = | = This Page gathers information about MRI processing (functional, using FS-FAST, and DTI, using FS and FSL tools) performed in the Language and Reading Research Lab (LRRL) = |
This Page gathers information about MRI processing (functional, using FS-FAST, and DTI, using FS and FSL tools) performed in the Language and Reading Research Lab (LRRL)
Please note that these commands are provided without guarantee and without support - Use at your own risks. Subject IDs are there as examples and should not be used (not that you could anyway). You should use these commands as guidelines, not as a script. They were provided as a response to many questions we have been receiving about the way we do things. It is CRITICAL that you know what each command and option does, please run the command names without any argument to check usage and meaning of options.
fMRI
#Move to analysis directory
cd tango<7/9>/<controls_analysis/ASD_analysis>
#set up env
nmrenv
setenv SUBJECTS_DIR /space/iggy/4/users/mody/recon_modygrp2/cherif
#Unpacking the data:
unpacksdcmdir -src <> -targ <> -scanonly file #to ID the run numbers.
--> results are saved in <subj>/unpack.log
unpacksdcmdir -src <> -targ <> -fsfast -run <04> bold NII f.nii -run <05> bold NII f.nii <add more runs as needed> -unpackerr
#Pre-process the data:
preproc-sess -s <subjid> -smout fmcsm5 -fwhm 5
Check Motion:
plot-twf-sess -s <> -mc
Check spiking:
spikedet-sess -s <>
--> exclude run or slices in question (use tpexclude)
cat ><subjID>/subjectname
<subjid>
<ctrl+d>
#Register functional to structural space (need the T1 to be recon'd):
fslregister-sess -s <subjid>
Check Registration:
tkregister-sess -s <subjid>
#Prepare model:
put together PAR files for subject in excel/Calc using the Presentation outputs and copy these to <subj>/bold/
cp <mody###_par_run1.txt><run#1>reason_par.par
and similarly for all par files in each bold/run directory
# Make Analysis
mkanalysis-sess -a pict-reas-dev-fir-noautostmdr-ter2 -TR 2.0 -paradigm reason_par.par -designtype event-related -funcstem fmcsm5 -nconditions 4 -timewindow 40 -prestim 4 -inorm -mcextreg -polyfit 2 -taumax 30 -TER 2 -noautostimdur -force
mkanalysis-sess -a pic-reas-dev-gf-autostim-ter2 -TR 2.0 -paradigm reason_par.par -designtype event-related -funcstem fmcsm5 -nconditions 4 -gammafit 2.25 1.25 -inorm -mcextreg -polyfit 2 -taumax 30 -TER 2 -autostimdur
#Make contrasts
mkcontrast-sess -analysis pic-reas-dev-fir-noautostmdr-ter2 -contrast omnibus -a 1 -a 2 -a 3 -c 0 -nosumconds
mkcontrast-sess -analysis pic-reas-dev-fir-noautostmdr-ter2 -contrast vs-v-fix -a 1 -c 0
mkcontrast-sess -analysis pic-reas-dev-fir-noautostmdr-ter2 -contrast perc-v-fix -a 2 -c 0
mkcontrast-sess -analysis pic-reas-dev-fir-noautostmdr-ter2 -contrast conc-v-fix -a 3 -c 0
mkcontrast-sess -analysis pic-reas-dev-fir-noautostmdr-ter2 -contrast vs-v-perc -a 1 -c 2
mkcontrast-sess -analysis pic-reas-dev-fir-noautostmdr-ter2 -contrast vs-v-conc -a 1 -c 3
mkcontrast-sess -analysis pic-reas-dev-fir-noautostmdr-ter2 -contrast perc-v-conc -a 2 -c 3
+same with rmprestim
mkcontrast-sess -analysis pic-reas-dev-gf-autostim-ter2 -contrast omnibus -a 1 -a 2 -a 3 -c 0 -nosumconds
mkcontrast-sess -analysis pic-reas-dev-gf-autostim-ter2 -contrast vs-v-fix -a 1 -c 0
mkcontrast-sess -analysis pic-reas-dev-gf-autostim-ter2 -contrast perc-v-fix -a 2 -c 0
mkcontrast-sess -analysis pic-reas-dev-gf-autostim-ter2 -contrast conc-v-fix -a 3 -c 0
mkcontrast-sess -analysis pic-reas-dev-gf-autostim-ter2 -contrast vs-v-perc -a 1 -c 2
mkcontrast-sess -analysis pic-reas-dev-gf-autostim-ter2 -contrast vs-v-conc -a 1 -c 3
mkcontrast-sess -analysis pic-reas-dev-gf-autostim-ter2 -contrast perc-v-conc -a 2 -c 3
mkanalysis-sess -a pic_reas_gf_nosmooth -TR 2.0 -paradigm reason _par.par -designtype event-related -funcstem fmc -nconditions 4 -gammafit 2.25 1 .25 -inorm -mcextreg -polyfit 2 -taumax 30 -TER 2 -autostimdur
mkanalysis-sess -a pic_reas_fir_nosmooth -TR 2.0 -paradigm reaso n_par.par -designtype event-related -funcstem fmc -nconditions 4 -timewindow 40 -prestim 4 -inorm -mcextreg -polyfit 2 -taumax 30 -TER 2 -noautostimdur
mkcontrast-sess -a pic_reas_gf_nosmooth -contrast omnibus -a 1 - a 2 -a 3 -c 0 -nosumconds
mkcontrast-sess -analysis pic_reas_gf_nosmooth -contrast omnibus -a 1 -a 2 -a 3 -c 0 -nosumconds
mkcontrast-sess -analysis pic_reas_gf_nosmooth -contrast vs-v-fi x -a 1 -c 0
mkcontrast-sess -analysis pic_reas_gf_nosmooth -contrast perc-v- fix -a 2 -c 0
mkcontrast-sess -analysis pic_reas_gf_nosmooth -contrast conc-v- fix -a 3 -c 0
mkcontrast-sess -analysis pic_reas_gf_nosmooth -contrast vs-v-pe rc -a 1 -c 2
mkcontrast-sess -analysis pic_reas_gf_nosmooth -contrast vs-v-co nc -a 1 -c 3
mkcontrast-sess -analysis pic_reas_gf_nosmooth -contrast perc-v- conc -a 2 -c 3
mkcontrast-sess -analysis pic_reas_fir_nosmooth -contrast omnibu s -a 1 -a 2 -a 3 -c 0 -nosumconds -rmprestim
mkcontrast-sess -analysis pic_reas_fir_nosmooth -contrast vs-v-f ix -a 1 -c 0 -rmprestim
mkcontrast-sess -analysis pic_reas_fir_nosmooth -contrast perc-v -fix -a 2 -c 0 -rmprestim
mkcontrast-sess -analysis pic_reas_fir_nosmooth -contrast conc-v -fix -a 3 -c 0 -rmprestim
mkcontrast-sess -analysis pic_reas_fir_nosmooth -contrast vs-v-p erc -a 1 -c 2 -rmprestim
mkcontrast-sess -analysis pic_reas_fir_nosmooth -contrast vs-v-c onc -a 1 -c 3 -rmprestim
mkcontrast-sess -analysis pic_reas_fir_nosmooth -contrast perc-v -conc -a 2 -c 3 -rmprestim
# Running the stats:
selxavg3-sess -a <analysis> -s <subj>
for gf, gf(nosmoothing), fir
or in seychelles:
foreach subject (mody028 mody050-new mody065 mody083 mody096 mody100 mody102 mody125 mody169 mody170 mody177 mody179 mody181)
foreach? pbsubmit -c "selxavg3-sess -a pic_reas_fir_nosmooth -s $subject" -m cherif
foreach? pbsubmit -c "selxavg3-sess -a pic_reas_gf_nosmooth -s $subject" -m cherif
foreach? end
#Viewing and saving statistical maps:
in each <subj>/bold/
mri_concat pic-reas-dev-gf-autostim-ter2/vs-v-fix/sig.nii pic-reas-dev-gf-autostim-ter2/perc-v-fix/sig.nii pic-reas-dev-gf-autostim-ter2/conc-v-fix/sig.nii pic-reas-dev-gf-autostim-ter2/vs-v-perc/sig.nii pic-reas-dev-gf-autostim-ter2/vs-v-conc/sig.nii pic-reas-dev-gf-autostim-ter2/perc-v-conc/sig.nii --o all.sig.nii
mkdir pics
cd pics
tksurfer <subj><hemi> inflated -overlay all.sig.nii
then in a separate terminal:
import <subj>_<hemi>_<contrast>_<view>.jpeg
and click on window.
#GROUP STATS:
Concatenate all subjects for a given analysis:
isxconcat-sess -s <subjid1> -s <subjid2>... -a <analysis> -all-contrasts -o <outputname> -hemis
cd <groupanalysis>/<analysis>/
FULL COHORT MEAN:
foreach cond (omnibus vs-v-fix perc-v-fix conc-v-fix vs-v-perc vs-v-conc perc-v-conc)
mri_glmfit --y $cond/lh.ces.nii --osgm --mask lh.mask.nii --glmdir $cond/lh.rfx.osgm --fwhm 5 --nii --surf fsaverage lh #for random effects ## do same for rh
mri_glmfit --y $cond/lh.ces.nii --osgm --mask lh.mask.nii --glmdir $cond/lh.wrfx.osgm --fwhm 5 --wls $cond/lh.cesvar.nii --nii --surf fsaverage lh #for weighted
random effects ## do same for rh
mri_glmfit --y $cond/lh.ces.nii --osgm --mask lh.mask.nii --glmdir $cond/lh.ffx.osgm --fwhm 5 --yffxvar $cond/lh.cesvar.nii --ffxdofdat ffxdot.dat --nii --surf fsaverage lh #for fixed effects ## do same for rh
mri_concat $cond/lh.rfx.osgm/osgm/sig.nii $cond/lh.wrfx.osgm/osgm/sig.nii $cond/lh.ffx.osgm/osgm/sig.nii -o $cond/lh.all.sig.nii
viewing results:
tksurfer fsaverage lh inflated -overlay <cond>/lh.all.sig.nii
#BETWEEN-GROUP STATS: ********************
isxconcat-sess -s <subjid1> -s <subjid2>... -a <analysis> -all-contrasts -o <outputname> -hemis
cd <groupdir>/<groupAnalysis>/<lowlevelanalysisname>
#create fsgdf as:
>>>>>>>
Title 4 kids vs. 4 adults
Class kids
Class adults
variables #can add variables such as age, RT etc.
input mody028 kids
input mody050 kids
input mody096 kids
input mody169 kids
input mody014 adults
input mody068 adults
input mody086 adults
input mody171 adults
>>>>>>>
# Make contrasts vectors
cat > kids_vs_adults.mat
1 -1
<ctrl+d>
cat > adults_vs_kids.mat
-1 1
<ctrl+d>
#Run stats:
foreach cond (vs-v-fix perc-v-fix conc-v-fix vs-v-perc vs-v-conc perc-v-conc)
mri_glmfit --y $cond/lh.ces.nii --fsgd <fsgdf> doss --C kids_vs_adults.mat --C adults_vs_kids.mat --mask lh.mask.nii --glmdir $cond/lh.rfx.osgm --nii --surf fsaverage lh #for random effects ## do same for rh
mri_glmfit --y $cond/lh.ces.nii --fsgd <fsgdf> doss --C kids_vs_adults.mat --C adults_vs_kids.mat --mask lh.mask.nii --glmdir $cond/lh.wrfx.osgm --wls $cond/lh.cesvar.nii --nii --surf fsaverage lh #for weighted random effects ## do same for rh
mri_glmfit --y $cond/lh.ces.nii --fsgd <fsgdf> doss --C kids_vs_adults.mat --C adults_vs_kids.mat --mask lh.mask.nii --glmdir $cond/lh.ffx.osgm --yffxvar $cond/lh.cesvar.nii --ffxdofdat ffxdot.dat --nii --surf fsaverage lh #for fixed effects ## do same for rh
CORRECTION FOR MULTIPLE COMPARISONS
random effects group map calculation:
pbsubmit -c "mri_glmfit --y $cond/lh.ces.nii --osgm --mask lh.mask.nii --glmdir $cond/lh.rfx.osgm --nii --fwhm 5 --surf fsaverage lh" -m cherif -o '-l nodes=1:opteron'
pbsubmit -c "mri_glmfit --y $cond/rh.ces.nii --osgm --mask rh.mask.nii --glmdir $cond/rh.rfx.osgm --nii --fwhm 5 --surf fsaverage rh" -m cherif -o '-l nodes=1:opteron'
then the same thing with a simulation:
pbsubmit -c "mri_glmfit --y $cond/lh.ces.nii --osgm --mask lh.mask.nii --glmdir $cond/lh.rfx.osgm --nii --fwhm 5 --surf fsaverage lh --sim perm 10000 1.3 $cond/lh_rfx_10000perm13" -m cherif -o '-l nodes=1:opteron'
pbsubmit -c "mri_glmfit --y $cond/rh.ces.nii --osgm --mask rh.mask.nii --glmdir $cond/rh.rfx.osgm --nii --fwhm 5 --surf fsaverage rh --sim perm 10000 1.3 $cond/rh_rfx_10000perm13" -m cherif -o '-l nodes=1:opteron'
Clustering: apply the csd in mri_durfcluster:
foreach cond ( vs-v-fix perc-v-fix conc-v-fix vs-v-perc vs-v-conc perc-v-conc )
foreach contrast (CTRL_grp HFA_grp CTRL-HFA HFA-CTRL)
mri_surfcluster --in $cond/lh.rfx/$contrast/sig.nii --csd $cond/lh_rfx_10000perm1.3-$contrast.csd --sum $cond/lh_rfx_10000csd_${contrast}-corrected.txt --o $cond/lh_rfx_10000csd_${contrast}-corrected.nii --cwsig $cond/lh_rfx_10000cwsig_${contrast}-corrected.nii
mri_surfcluster --in $cond/rh.rfx/$contrast/sig.nii --csd $cond/rh_rfx_10000perm13-$contrast.csd --sum $cond/rh_rfx_10000csd_${contrast}-corrected.txt --o $cond/rh_rfx_10000csd_${contrast}-corrected.nii --cwsig $cond/rh_rfx_10000cwsig_${contrast}-corrected.nii
end
If you want the corrected p values then use the --cwsig output. Each voxel in a cluster will have the same p-value (the p-value of the cluster).
Can do the same with Weighted random effects, or fixed effect...
foreach cond ( vs-v-fix perc-v-fix conc-v-fix vs-v-perc vs-v-conc perc-v-conc )
pbsubmit -c "mri_glmfit --y $cond/lh.ces.nii --osgm --mask lh.mask.nii --glmdir $cond/lh.wrfx.osgm --wls $cond/lh.cesvar.nii --nii --fwhm 5 --surf fsaverage lh" -m cherif -o '-l nodes=1:opteron'
pbsubmit -c "mri_glmfit --y $cond/rh.ces.nii --osgm --mask rh.mask.nii --glmdir $cond/rh.wrfx.osgm --wls $cond/lh.cesvar.nii --nii --fwhm 5 --surf fsaverage rh" -m cherif -o '-l nodes=1:opteron'
pbsubmit -c "mri_glmfit --y $cond/lh.ces.nii --osgm --mask lh.mask.nii --glmdir $cond/lh.wrfx.osgm --wls $cond/lh.cesvar.nii --nii --surf fsaverage lh --sim perm 2000 1.3 $cond/lh_wrfx_2000perm13" -m cherif -o '-l nodes=1:opteron'
pbsubmit -c "mri_glmfit --y $cond/rh.ces.nii --osgm --mask rh.mask.nii --glmdir $cond/rh.wrfx.osgm --wls $cond/lh.cesvar.nii --nii --fwhm 5 --surf fsaverage rh --sim perm 2000 1.3 $cond/rh_wrfx_2000perm13" -m cherif -o '-l nodes=1:opteron'
viewing results:
mri_concat vs-v-fix/lh_rfx_10000csd_corrected.nii etc.
tksurfer fsaverage lh inflated -overlay <cond>/concat_all_sigcorr.nii
# Practical Example:
Group stats 10CTRL_vs_10HFA
foreach cond ( vs-v-fix perc-v-fix conc-v-fix vs-v-perc vs-v-conc perc-v-conc )
foreach? mri_glmfit --y $cond/lh.ces.nii --fsgd 10CTRL-10HFA_fsgdf_nocov doss --C HFA-only.mat --C CTRL-only.mat --C CTRL_vs_HFA_nocov.mat --mask lh.mask.nii --glmdir $cond/lh.rfx.nocov --fwhm 5 --nii --surf fsaverage lh
foreach? mri_glmfit --y $cond/rh.ces.nii --fsgd 10CTRL-10HFA_fsgdf_nocov doss --C HFA-only.mat --C CTRL-only.mat --C CTRL_vs_HFA_nocov.mat --mask rh.mask.nii --glmdir $cond/rh.rfx.nocov --fwhm 5 --nii --surf fsaverage rh
foreach? end
foreach cond ( vs-v-fix perc-v-fix conc-v-fix vs-v-perc vs-v-conc perc-v-conc )
mri_glmfit --y $cond/lh.ces.nii --fsgd 10CTRL-10HFA_fsgdf_nocov doss --C HFA-only.mat --C CTRL-only.mat --C CTRL_vs_HFA_nocov.mat --mask lh.mask.nii --glmdir $cond/lh.wrfx.nocov --wls $cond/lh.cesvar.nii --nii --fwhm 5 --surf fsaverage lh
mri_glmfit --y $cond/rh.ces.nii --fsgd 10CTRL-10HFA_fsgdf_nocov doss --C HFA-only.mat --C CTRL-only.mat --C CTRL_vs_HFA_nocov.mat --mask rh.mask.nii --glmdir $cond/rh.wrfx.nocov --wls $cond/lh.cesvar.nii --nii --fwhm 5 --surf fsaverage rh
pbsubmit -c "mri_glmfit --y $cond/lh.ces.nii --fsgd 10CTRL-10HFA_fsgdf_nocov doss --C HFA-only.mat --C CTRL-only.mat --C CTRL_vs_HFA_nocov.mat --mask lh.mask.nii --glmdir $cond/lh.wrfx.nocov --wls $cond/lh.cesvar.nii --nii --surf fsaverage lh --sim perm 2000 1.3 $cond/lh_wrfx_2000perm13" -m cherif -o '-l nodes=1:opteron'
pbsubmit -c "mri_glmfit --y $cond/rh.ces.nii --fsgd 10CTRL-10HFA_fsgdf_nocov doss --C HFA-only.mat --C CTRL-only.mat --C CTRL_vs_HFA_nocov.mat --mask rh.mask.nii --glmdir $cond/rh.wrfx.nocov --wls $cond/lh.cesvar.nii --nii --fwhm 5 --surf fsaverage rh --sim perm 2000 1.3 $cond/rh_wrfx_2000perm13" -m cherif -o '-l nodes=1:opteron'
mri_glmfit --y $cond/lh.ces.nii --fsgd 10CTRL-10HFA_fsgdf_nocov doss --C HFA-only.mat --C CTRL-only.mat --C CTRL_vs_HFA_nocov.mat --mask lh.mask.nii --glmdir $cond/lh.ffx.nocov --fwhm 5 --yffxvar $cond/lh.cesvar.nii --ffxdofdat ffxdof.dat --nii --surf fsaverage lh
mri_glmfit --y $cond/rh.ces.nii --fsgd 10CTRL-10HFA_fsgdf_nocov doss --C HFA-only.mat --C CTRL-only.mat --C CTRL_vs_HFA_nocov.mat --mask rh.mask.nii --glmdir $cond/rh.ffx.nocov --fwhm 5 --yffxvar $cond/rh.cesvar.nii --ffxdofdat ffxdof.dat --nii --surf fsaverage rh
#OMNIBUS CALCULATION AND ROI DRAWING:
mri_glmfit --y omnibus/lh.ces.000.nii --osgm --mask lh.mask.nii --glmdir omnibus/lh.ces.000.rfx.osgm --nii --surf fsaverage lh --fwhm 5
mri_glmfit --y omnibus/lh.ces.001.nii --osgm --mask lh.mask.nii --glmdir omnibus/lh.ces.001.rfx.osgm --nii --surf fsaverage lh --fwhm 5
mri_glmfit --y omnibus/lh.ces.002.nii --osgm --mask lh.mask.nii --glmdir omnibus/lh.ces.002.rfx.osgm --nii --surf fsaverage lh --fwhm 5
mri_glmfit --y omnibus/rh.ces.000.nii --osgm --mask rh.mask.nii --glmdir omnibus/rh.ces.000.rfx.osgm --nii --surf fsaverage rh --fwhm 5
mri_glmfit --y omnibus/rh.ces.001.nii --osgm --mask rh.mask.nii --glmdir omnibus/rh.ces.001.rfx.osgm --nii --surf fsaverage rh --fwhm 5
mri_glmfit --y omnibus/rh.ces.002.nii --osgm --mask rh.mask.nii --glmdir omnibus/rh.ces.002.rfx.osgm --nii --surf fsaverage rh --fwhm 5
Each 000, 001, 002 was deconstructed into individual t-tests from the f-test omnibus (i.e. each one is now an EV/condition)
to put them back together, take the best p-value of each test at each voxel and bonferroni correct (i.e multiply by 3)
In practice:
mri_concat sig1.nii sig2.nii sig3.nii --max-bonfcor --o sig.max.bonfcor.nii --abs
-->
mri_concat lh.ces.000.rfx.osgm/osgm/sig.nii lh.ces.001.rfx.osgm/osgm/sig.nii lh.ces.002.rfx.osgm/osgm/sig.nii --max-bonfcor --o lh.rfx.sigmax.bonfcor.nii --abs
mri_concat rh.ces.000.rfx.osgm/osgm/sig.nii rh.ces.001.rfx.osgm/osgm/sig.nii rh.ces.002.rfx.osgm/osgm/sig.nii --max-bonfcor --o rh.rfx.sigmax.bonfcor.nii --abs
#DRAWING ROIs:
right click to reset points
click around the pic to make boundary points
click make closed boundary
click in the middle of it and then click custom fill --> up to and including paths
save label as...
repeat after right clicking brain.
all in 13K_fMRI_gf/ROIs/
#Extracting timecourses for all subjects in predefined ROIs:
cp 13K_fMRI_gf/ROIs/* .
forach label (*.label)
func2roi-sess -roidef roi_$label -analysis pic-reas-dev-fir-noautostmdr-ter2 -labelfile $label -labelspace tal -s mody028 -s mody050-new -s mody065 -s mody083 -s mody096 -s mody100 -s mody102 -s mody125 -s mody169 -s mody170 -s mody177 -s mody179 -s mody181
end
forach label (*.label)
func2roi-sess -roidef roi_$label -analysis pic_reas_fir_nosmooth -labelfile $label -labelspace tal -s mody028 -s mody050-new -s mody065 -s mody083 -s mody096 -s mody100 -s mody102 -s mody125 -s mody169 -s mody170 -s mody177 -s mody179 -s mody181
end
foreach label ( rh*.label )
foreach? pbsubmit -c "func2roi-sess -roidef roi_$label -analysis pic-reas-dev-fir-noautostmdr-ter2 -labelfile $label -labelspace tal -s mody028 -s mody050-new -s mody065 -s mody083 -s mody096 -s mody100" foreach? pbsubmit -c "func2roi-sess -roidef roi_$label -analysis pic-reas-dev-fir-noautostmdr-ter2 -labelfile $label -labelspace tal -s mody102 -s mody125 -s mody169 -s mody170 -s mody177 -s mody179 -s mody181" foreach? pbsubmit -c "func2roi-sess -roidef roi_$label -analysis pic-reas-dev-fir-noautostmdr-ter2 -labelfile $label -labelspace tal -s mody175 -s mody187 -s mody188 -s mody189 -s mody191 -s mody192" foreach? pbsubmit -c "func2roi-sess -roidef roi_$label -analysis pic-reas-dev-fir-noautostmdr-ter2 -labelfile $label -labelspace tal -s mody195 -s mody197 -s mody198 -s mody199"
mkdir roisum
foreach label (*.label)
roisummary-sess -sumfile roisum/$label -roidef roi_$label -analysis pic-reas-dev-fir-noautostmdr-ter2 -s mody028 -s mody050-new -s mody065 -s mody083 -s mody096 -s mody100 -s mody102 -s mody125 -s mody169 -s mody170 -s mody177 -s mody179 -s mody181
end
foreach label (*.label)
roisummary-sess -sumfile roisum_nosmooth/$label -roidef roi_$label -analysis pic_reas_fir_nosmooth -s mody028 -s mody050-new -s mody065 -s mody083 -s mody096 -s mody100 -s mody102 -s mody125 -s mody169 -s mody170 -s mody177 -s mody179 -s mody181
end
this gives a text file as described in the wiki for each label
from there, line 9 is the number of conditions, line 10 the number of estimates per condition, and the following rows are the number estimates. each of these estimates is divided by the offset of that column (line4) and multiplied by 100 to ge the %signal change.
or on Seychelles:
foreach label ( lh*.label )
foreach? pbsubmit -c "roisummary-sess -sumfile roisum/CTRL_lh_$label -roidef roi_$label -analysis pic-reas-dev-fir-noautostmdr-ter2 -s mody050-new -s mody065 -s mody083 -s mody096 -s mody100 -s mody125 -s mody169 -s mody170 -s mody177 -s mody179"
foreach? end
foreach label ( lh*.label )
foreach? pbsubmit -c "roisummary-sess -sumfile roisum/HFA_lh_$label -roidef roi_$label -analysis pic-reas-dev-fir-noautostmdr-ter2 -s mody175 -s mody187 -s mody188 -s mody189 -s mody191 -s mody192 -s mody195 -s mody197 -s mody198 -s mody199"
foreach? end
DTI
unpacksdcmdir -src /space/archive/95/siemens/<TrioTim-35006-20060923-084459-893000> -targ /space/oscar/6/users/DTI/<RPM_001> -run <11> dti NII <rpm-001>-dti.nii -unpackerr
cd <subjid>/dti/
fslswapdim <run#>/<in> x -y z <out>
fslorient -forceradiological <out>
where in is <subj>_dti.nii and out is swap_<subj>_dti.nii
# EDDY CURRENT CORRECTION:
#GUI
>>Fdt
GUI: input volume swap<rpm-001>_dti.hdr, output volume data, ref volume 0
#LINE COMMAND:
>> eddy_correct <input_volume><output_name> 0
#where 0 is the reference volume for registration, this needs to be a volume where no DTI gradient was applied (in our case one of the first 10)
#MAKE SURE TO USE "DATA" AS THE OUTPUT FILENAME, THSI SIMPLIFIES LATER STEPS BY PROVIDING AN EXPECTED STRUCTURE
#CONSTRUCTING THE BRAIN MASK:
# FIRST ISOLATE THE VOLUMES:
fslsplit swap<rpm-001-dti>
mv vol0000.nii.gz nodif-vol0000.nii.gz
rm vol*
bet nodif-vol000o.nii.gz nodif_brain.nii.gz -m -f 0.3
fslview nodif_brain #check brain extraction
# so we get nodif_brain and nodif_brain_mask
#ENSURE YOU HAVE the bvals and bvecs files in the appropriate format - This should not change for a given acquisition, so we're all set with our master files.
>>cp <path>/bvals .
>>cp <path>/bvecs .
#HERE CHECK THE DIRECTORY STRUCTURE IS CORRECT:
>> bedpostx_datacheck .
#CALCULATION OF DIFFUSION TENSOR MODEL AT EACH VOXEL:
>>dtifit -k data.nii.gz -o dtifit -m nodif_brain_mask.nii.gz -r bvecs -b bvals
# if you have kept the filenames suggested, then the above line can be used AS SUCH
#BUILDING UP DISTRIBUTIONS ON DIFFUSION PARAMETERS AT EACH VOXEL - YOU DO NOT NEED TO HAVE RUN DTIFIT FOR THIS - THIS WILL CREATE THE FILES NECESSARY FOR TRACTOGRAPHY - ANASTASIA'S NEW BEDPOSTX code:
#the input directory option can only be used if the filename structure above was used, otherwise, need to enter all the inputs separately.
[seychelles:012] (nmr-dev-env) bedpostx_seychelles .
# this will submit the jobs for me, i don't need to call pbsubmit
Type /homes/6/cherif/iggy2/controls_analysis/mody102/dti/012.bedpostX/monitor on seychelles to show progress.
Type /homes/6/cherif/iggy2/controls_analysis/mody102/dti/012.bedpostX/cputime on seychelles to show CPU time of completed jobs.
Type /homes/6/cherif/iggy2/controls_analysis/mody102/dti/012.bedpostX/cancel on seychelles to terminate all the queued tasks.
You will get an email at the end of the post-processing stage.
#TRACKING
#transform labels into binary ROIs:
Draw ROIs on omnibus and save in
20_omnibus_labels
then
#Morph all labels to individual surface space
foreach label (lh*.label)
foreach? foreach subject (mody050-new mody065 mody083 mody096 mody100 mody125 mody169 mody170 mody177 mody179 mody175 mody187 mody188 mody189 mody191 mody192 mody195 mody197 mody198 mody199)
foreach? mri_label2label --srclabel $label --srcsubject fsaverage --trgsubject $subject --trglabel $label --regmethod surface --hemi lh
foreach? end
foreach? end
foreach label (rh*.label)
foreach subject (mody028 mody050-new mody065 mody083 mody096 mody100 mody102 mody125 mody169 mody170 mody177 mody179 mody181)
mri_label2label --srclabel $label --srcsubject fsaverage --trgsubject $subject --trglabel $label --regmethod surface --hemi rh
end
end
foreach subject (mody*)
mv $subject/dti/T1_ROIs $subject/dti/T1_ROIs_old
mkdir $subject/dti/T1_ROIs
end
#Transform surface labels into volume ROIs
foreach label (lh.IPS.label lh_supFrontalS.label lh_supPrecentralS.label lh.TPJ.label)
foreach? foreach subject (mody*)
foreach? mri_label2vol --label $subject/label/$label --temp $subject/mri/brain.mgz --fillthresh 0.5 --o /space/tango/9/users/mody/cherif/HFA/$subject/dti/T1_ROIs/$label.nii
foreach? end
#BUT FOR MORE ACCURATE TRACKING, USE LARGER RIBBON FOR SEED MASKS:
#This uses a ribbon 2.5mm into the white matter space
foreach label ( lh.fusiform.label lh.post-cingulate.label lh_supPrecentralS.label lh.IFS.label lh.STS.label lh.supramarginal.label lh.IPS.label lh_supFrontalS.label lh.TPJ.label )
foreach? foreach subject ( mody050-new mody065 mody083 mody096 mody100 mody125 mody169 mody170 mody171 mody175 mody177 mody179 mody187 mody188 mody189 mody191 mody192 mody195 mody197 mody198 mody199)
foreach? mri_label2vol --label $SUBJECTS_DIR/$subject/label/$label --temp $SUBJECTS_DIR/$subject/mri/brain.mgz --proj abs -2.5 0 .1 --o $subject/dti/T1_ROIs/$label.proj2.5.nii --subject $subject --hemi lh
foreach? end
foreach? end
#Need the T1_brain volume as a nifti file:
foreach subject (mody191 mody192 mody195)
mri_convert $SUBJECTS_DIR/$subject/mri/brain.mgz $subject/dti/T1_brain.nii
end
#Calculate resampling matrix from DTI space to T1 space and invert it to get T12FA to be able to transform all volume ROIs to FA space for tracking
foreach subject (mody*)
pbsubmit -c "flirt -in $subject/dti/dtifit_FA -ref $subject/dti/T1_brain -omat $subject/dti/FA2T1brain.mat -out $subject/dti/FA2T1brain.nii.gz -searchrx -180 180 -searchry -180 180 -searchrz -180 180"
end
now back out of seychelles:
foreach subject (mody*)
convert_xfm -omat $subject/T12FA.mat -inverse $subject/FA2T1.mat
end
and into seychelles again:
#Apply transform to all subjects for each label
foreach fROI ( lh.fusiform.label.proj2.5.nii lh_supFrontalS.label.proj2.5.nii rh.fusiform.label.proj2.5.nii lh.IFS.label.proj2.5.nii lh_supPrecentralS.label.proj2.5.nii rh.IPS.label.proj2.5.nii lh.IPS.label.proj2.5.nii lh.supramarginal.label.proj2.5.nii rh.post-cingulate.label.proj2.5.nii lh.post-cingulate.label.proj2.5.nii lh.TPJ.label.proj2.5.nii rh.supramarginal.label.proj2.5.nii lh.STS.label.proj2.5.nii rh.ant-cingulate.label.proj2.5.nii rh.TPJ.label.proj2.5.nii )
foreach? foreach subject (mody050-new mody065 mody083 mody096 mody100 mody125 mody169 mody170 mody171 mody175 mody177 mody179 mody187 mody188 mody189 mody191 mody192 mody195 mody197 mody198 mody199)
foreach? pbsubmit -c "flirt -in $subject/dti/T1_ROIs/$fROI -applyxfm -init $subject/dti/T1brain2FA.mat -out $subject/dti/DTI_ROIs/$fROI -ref $subject/dti/dtifit_FA.nii.gz"
foreach? end
foreach? end
foreach fROI ( lh.fusiform.label.proj2.5.nii lh_supFrontalS.label.proj2.5.nii rh.fusiform.label.proj2.5.nii lh.IFS.label.proj2.5.nii lh_supPrecentralS.label.proj2.5.nii rh.IPS.label.proj2.5.nii lh.IPS.label.proj2.5.nii lh.supramarginal.label.proj2.5.nii rh.post-cingulate.label.proj2.5.nii lh.post-cingulate.label.proj2.5.nii lh.TPJ.label.proj2.5.nii rh.supramarginal.label.proj2.5.nii lh.STS.label.proj2.5.nii rh.ant-cingulate.label.proj2.5.nii rh.TPJ.label.proj2.5.nii )
foreach? foreach subject (mody050-new mody065 mody083 mody096 mody100 mody125 mody169 mody170 mody175 mody177 mody179 mody187 mody188 mody189 mody191 mody192 mody195 mody197 mody198 mody199)
foreach? fslmaths $subject/dti/DTI_ROIs/$fROI -bin $subject/dti/DTI_ROIs/$fROI
foreach? end
foreach? end
#TO RUN PROBTRACKX IN STD MNI SPACE:
foreach label
mri_label2vol --label $label --temp $SUBJECTS_DIR/fsaverage/mri/brain.mgz --proj abs -2.5 0 .1 --o ${label}_vol.nii --subject fsaverage --hemi <l/r>h
Need std2DTI transforms as a loop for all subjects:
flirt FA to T1:
as above
flirt T1 to fsaverage:
foreach subject (mody*)
foreach? flirt -in $subject/dti/T1_brain.nii.gz -ref 20_omnibus_labels/fsaverage_brain.nii.gz -omat $subject/dti/T1brain2fsavg.mat -out $subject/dti/T1brain2fsavg.nii.gz -searchrx -180 180 -searchry -180 180 -searchrz -180 180
foreach? end
convert_xfm -omat FA2std -concat FA2T1 T12std
convert_xfm -omat std2Fa -inverse FA2std
#Make Exclusion Mask:
In fslview, on MNI_brain volume, draw an hemispheric mask, on the midline, masking all but the corpus callosum, then register this mask onto each subject using FLIRT with the MNI2FA.mat.
Open each mask on top of dtifit_new_FA and make sure it's in the right place, adding/deleting if necessary.
The recon already created lh.cortex.label and rh.cortex.label. Morph these onto the T1 volume, with -proj abs 3 3.5 .1 to keep it inside the grey matter when you morph onto FA maps.
Flirt the cortex maps onto FA space.
fslmaths to add the lh, rh cortices, and hemispheric masks.
fslmaths -bin to binarize the exclusion mask as exclusionmasktotal.nii.gz
#Run all probabilistic trackings between areas of interest:
- foreach subject ( mody028 mody050-new mody065 mody073_new mody096 mody100 mody125 mody126_new mody169 mody177 mody179 mody181 mody175 mody187 mody188 mody189 mody192 mody194 mody195 mody197 mody199 mody201 mody202 mody203 )
pbsubmit -c "probtrackx --mode=seedmask -x $subject/dti/DTI_ROIs/lh.IFS.label.proj2.5.nii --avoid=$subject/dti/DTI_ROIs/exclusionmask.nii.gz -c 0.01 -S 2000 --steplength=0.5 -P 5000 --stop=$subject/dti/DTI_ROIs/lh.STS.label.proj2.5.nii --forcedir --opd -s $subject/dti/dti.bedpostX/merged -m $subject/dti/dti.bedpostX/nodif_brain_mask --dir=$subject/dti/dti.bedpostX/DTItrax_lhIFS2lharcuate2lhSTS_new_noloopcheck --waypoints=$subject/dti/DTI_ROIs/waypts_arcuate-STS.txt"
pbsubmit -c "probtrackx --mode=seedmask -x $subject/dti/DTI_ROIs/lh.IPS.label.proj2.5.nii --avoid=$subject/dti/DTI_ROIs/exclusionmask.nii.gz -l -c 0.01 -S 2000 --steplength=0.5 -P 5000 --stop=$subject/dti/DTI_ROIs/lh.occipital.label.proj2.5.nii --forcedir --opd -s $subject/dti/dti.bedpostX/merged -m $subject/dti/dti.bedpostX/nodif_brain_mask --dir=$subject/dti/dti.bedpostX/DTItrax_lhIPS2lhOcc_new --waypoints=$subject/dti/DTI_ROIs/lh.occipital.label.proj2.5.nii"
pbsubmit -c "probtrackx --mode=seedmask -x $subject/dti/DTI_ROIs/lh.IPS.label.proj2.5.nii --avoid=$subject/dti/DTI_ROIs/exclusionmask.nii.gz -l -c 0.01 -S 2000 --steplength=0.5 -P 5000 --stop=$subject/dti/DTI_ROIs/rh.IPS.label.proj2.5.nii --forcedir --opd -s $subject/dti/dti.bedpostX/merged -m $subject/dti/dti.bedpostX/nodif_brain_mask --dir=$subject/dti/dti.bedpostX/DTItrax_lhIPS2rhIPS_new --waypoints=$subject/dti/DTI_ROIs/rh.IPS.label.proj2.5.nii"
pbsubmit -c "probtrackx --mode=seedmask -x $subject/dti/DTI_ROIs/rh.occipital.label.proj2.5.nii --avoid=$subject/dti/DTI_ROIs/exclusionmask.nii.gz -l -c 0.01 -S 2000 --steplength=0.5 -P 5000 --stop=$subject/dti/DTI_ROIs/rh.IPS.label.proj2.5.nii --forcedir --opd -s $subject/dti/dti.bedpostX/merged -m $subject/dti/dti.bedpostX/nodif_brain_mask --dir=$subject/dti/dti.bedpostX/DTItrax_rhOcc2rhIPS_new --waypoints=$subject/dti/DTI_ROIs/rh.IPS.label.proj2.5.nii"
pbsubmit -c "probtrackx --mode=seedmask -x $subject/dti/DTI_ROIs/CC.nii --avoid=$subject/dti/DTI_ROIs/exclusionmask.nii.gz -l -c 0.01 -S 2000 --steplength=0.5 -P 5000 --forcedir --opd -s $subject/dti/dti.bedpostX/merged -m $subject/dti/dti.bedpostX/nodif_brain_mask --dir=$subject/dti/dti.bedpostX/DTItrax_CC_new"
Normalizing pathways:
We have tracked from A to B and from B to A, now to get the probability that either direction pathway passes through voxel x, normalized to waytotals (and therefore also allowing comparing these probabilities across subjects), the bottom line is to fslmaths your way to:
p = [p( [A->x->B] )]/waytotal(a->b) + [p( [A<-x<-B] )]/waytotal(b->a) - ([p( [A->x->B] )]/waytotal(a->b)) * ([p( [A<-x<-B] )]/waytotal(b->a))
(assuming independance, which is clearly an approximation, but one that is hard to get around)
SO:
foreach trak (DTItrax_lhIPS2lhIFS_new DTItrax_lhIPS2lhFG_new DTItrax_lhIPS2lhSTS_new DTItrax_lhOcc2lhIFS_new DTItrax_lhIFS2lhIPS_new DTItrax_lhSTS2lhIPS_new DTItrax_lhIFS2lhOcc_new DTItrax_lhTPJ2lhIPS_new DTItrax_lhIPS2lhTPJ_new DTItrax_lhFG2lhIPS_new DTItrax_lhFG2lhIFS_new DTItrax_lhIFS2lhFG_new DTItrax_lhTPJ2lhFG_new DTItrax_lhFG2lhTPJ_new DTItrax_lhTPJ2lharcuate2lhIFS_new DTItrax_lhIFS2lharcuate2lhTPJ_new DTItrax_lhSTS2lharcuate2lhIFS_new DTItrax_lhIFS2lharcuate2lhSTS_new)
foreach subject ( mody028 mody050-new mody065 mody073_new mody096 mody100 mody125 mody126_new mody169 mody177 mody179 mody181 mody175 mody187 mody188 mody189 mody192 mody194 mody195 mody197 mody199 mody201 mody202 mody203 )
set w1=head -1 $subject/dti/dti.bedpostX/$trak/waytotal | cut -f1 -d''
fslmaths $subject/dti/dti.bedpostX/$trak/fdt_paths -div $w1 $subject/dti/dti.bedpostX/$trak/fdt_paths_condprob.nii.gz -odt float
end
end
Checking distribution of tracks:
foreach trak (DTItrax_lhIPS2lhIFS_new DTItrax_lhIPS2lhFG_new DTItrax_lhIPS2lhSTS_new DTItrax_lhOcc2lhIFS_new DTItrax_lhIFS2lhIPS_new DTItrax_lhSTS2lhIPS_new DTItrax_lhIFS2lhOcc_new DTItrax_lhTPJ2lhIPS_new DTItrax_lhIPS2lhTPJ_new DTItrax_lhFG2lhIPS_new DTItrax_lhFG2lhIFS_new DTItrax_lhIFS2lhFG_new DTItrax_lhTPJ2lhFG_new DTItrax_lhFG2lhTPJ_new DTItrax_lhTPJ2lharcuate2lhIFS_new DTItrax_lhIFS2lharcuate2lhTPJ_new DTItrax_lhSTS2lharcuate2lhIFS_new DTItrax_lhIFS2lharcuate2lhSTS_new)
echo $trak >> Trax/histograms.txt
foreach subject ( mody028 mody050-new mody065 mody073_new mody096 mody100 mody125 mody126_new mody169 mody177 mody179 mody181 mody175 mody187 mody188 mody189 mody192 mody194 mody195 mody197 mody199 mody201 mody202 mody203 )
echo $subject >> Trax/histograms.txt
fslstats $subject/dti/dti.bedpostX/$trak/fdt_paths_condprob.nii.gz -H 24 0 1.2 >> Trax/histograms.txt
end
end
#Now import $trak.txt to Calc/Excel and plot the histogram --> If the histogram does not have a 1/x type of shape, then this fdt_paths is rejected from future processing.
Once You decide which fdt_paths to drop, you can compute fdt_paths_OR, as one path, or the Unionof 2 paths, as below:
foreach trak1 (DTItrax_lhOcc2lhIFS_new)
foreach trak2 (DTItrax_lhIFS2lhOcc_new)
foreach subject ( mody028 mody050-new mody065 mody073_new mody096 mody100 mody125 mody126_new mody169 mody177 mody179 mody181 mody175 mody187 mody188 mody189 mody192 mody194 mody195 mody197 mody199 mody201 mody202 mody203 )
fslmaths $subject/dti/dti.bedpostX/$trak1/fdt_paths_condprob.nii.gz -mul $subject/dti/dti.bedpostX/$trak2/fdt_paths_condprob.nii.gz $subject/dti/dti.bedpostX/$trak1/fdt_paths_intersect.nii.gz -odt float
fslmaths $subject/dti/dti.bedpostX/$trak1/fdt_paths_condprob.nii.gz -add $subject/dti/dti.bedpostX/$trak2/fdt_paths_condprob.nii.gz -sub $subject/dti/dti.bedpostX/$trak1/fdt_paths_intersect.nii.gz $subject/dti/dti.bedpostX/$trak1/fdt_paths_Union.nii.gz -odt float
end
Now manually go in and replace fdt_paths_Union with the single direction files as necessary.
Then:
foreach trak1 (DTItrax_lhOcc2lhIFS_new)
foreach subject ( mody028 mody050-new mody065 mody073_new mody096 mody100 mody125 mody126_new mody169 mody177 mody179 mody181 mody175 mody187 mody188 mody189 mody192 mody194 mody195 mody197 mody199 mody201 mody202 mody203 )
fslmaths $subject/dti/dti.bedpostX/$trak1/fdt_paths_Union.nii.gz -thr 0.05 -bin $subject/dti/dti.bedpostX/$trak1/fdt_paths_Union_bin0.05.nii.gz end
mkdir Trax/lhFG2lhIFS_OR
mkdir Trax/lhFG2lhIFS_OR_bin0.05
Extracting values:
foreach trak1 (DTItrax_lhOcc2lhIFS_new)
foreach subject ( mody028 mody050-new mody065 mody073_new mody096 mody100 mody125 mody126_new mody169 mody177 mody179 mody181 mody175 mody187 mody188 mody189 mody192 mody194 mody195 mody197 mody199 mody201 mody202 mody203 )
fslstats $subject/dti/dti.bedpostX/$trak1/fdt_paths_Union.nii.gz -l 0.05 -M
end
echo FA
foreach subject ( mody028 mody050-new mody065 mody073_new mody096 mody100 mody125 mody126_new mody169 mody177 mody179 mody181 mody175 mody187 mody188 mody189 mody192 mody194 mody195 mody197 mody199 mody201 mody202 mody203 )
fslstats $subject/dti/dtifit_new_FA.nii.gz -k $subject/dti/dti.bedpostX/$trak1/fdt_paths_Union_bin0.05.nii.gz -M
end
end
Extracting probability-weighted FA means:
foreach trak1 (DTItrax_lhOcc2lhIFS_new)
foreach subject ( mody028 mody050-new mody065 mody073_new mody096 mody100 mody125 mody126_new mody169 mody177 mody179 mody181 mody175 mody187 mody188 mody189 mody192 mody194 mody195 mody197 mody199 mody201 mody202 mody203 )
fslmaths $subject/dti/dti.bedpostX/$trak1/fdt_paths_Union.nii.gz -thr 0.05 $subject/dti/dti.bedpostX/$trak1/fdt_paths_Union_0.05.nii.gz
fslmaths $subject/dti/dtifit_new_FA.nii.gz -mul $subject/dti/dti.bedpostX/$trak1/fdt_paths_Union_0.05.nii.gz $subject/dti/dti.bedpostX/$trak1/weighted_FA_0.05.nii.gz
fslstats $subject/dti/dti.bedpostX/$trak1/weighted_FA_0.05.nii.gz -M
end
Making group Maps of pathway in MNI space:
foreach subject ( mody028 mody050-new mody065 mody073_new mody096 mody100 mody125 mody126_new mody169 mody177 mody179 mody181 mody175 mody187 mody188 mody189 mody192 mody194 mody195 mody197 mody199 mody201 mody202 mody203 )
flirt -in $subject/dti/dti.bedpostX/DTItrax_lhOcc2lhIFS_new/fdt_paths_Union.nii.gz -applyxfm -init $subject/dti/FA2MNI.mat -ref MNI_brain.nii.gz -out Trax/lhOcc2lhIFS_Union/${subject}.nii.gz
fslmaths Trax/lhOcc2lhIFS_Union/${subject}.nii.gz -bin Trax/lhOcc2lhIFS_Union/${subject}.bin.nii.gz
end
foreach subject ( mody028 mody050-new mody065 mody073_new mody096 mody100 mody125 mody126_new mody169 mody177 mody179 mody181 mody175 mody187 mody188 mody189 mody192 mody194 mody195 mody197 mody199 mody201 mody202 mody203 )
flirt -in $subject/dti/dti.bedpostX/DTItrax_lhFG2lhIFS_new/fdt_paths_OR_bin0.05.nii.gz -applyxfm -init $subject/dti/FA2MNI.mat -ref MNI_brain.nii.gz -out Trax/lhFG2lhIFS_OR_bin0.05/${subject}.nii.gz
fslmaths Trax/lhFG2lhIFS_OR_bin0.05/${subject}.nii.gz -bin Trax/lhFG2lhIFS_OR_bin0.05/${subject}.nii.gz
end
cd Trax/lhIPS2lhIFS_OR_bin0.05
fslmaths mody028.nii.gz -add mody050-new.nii.gz -add mody065.nii.gz -add mody073_new.nii.gz -add mody096.nii.gz -add mody100.nii.gz -add mody125.nii.gz -add mody126_new.nii.gz -add mody169.nii.gz -add mody177.nii.gz -add mody179.nii.gz -add mody181.nii.gz CTRLS.nii.gz
fslmaths mody175.nii.gz -add mody187.nii.gz -add mody188.nii.gz -add mody189.nii.gz -add mody192.nii.gz -add mody194.nii.gz -add mody195.nii.gz -add mody197.nii.gz -add mody199.nii.gz -add mody201.nii.gz -add mody202.nii.gz -add mody203.nii.gz HFA.nii.gz
#Extracting values from the pathway distributions obtained
Now need to normalize the tracts so they make sense for group comparison
things to take into accout:
- number of voxels in each seed
- total number of streamlines (waytotal)
- threshold tract based on %(waytotal)
total number of tracts sent out = #samples (5000) x #voxels in seed (=fslstats -V) = waytotal if no rejection, or >waytotal if paths were rejected
NORMALIZING TO WAYTOTAL: fdt_paths/waytotal for each subject gives equivalent voxels intensity
so:
#Create a file with all waytotal values for each subject
foreach s (mody*)
foreach? more $s/dti/dti.bedpostX//probtrackx_lhTPJ2lhsupPrecentral/waytotal >><track>_<groups>_waytotal.txt
foreach? end
we want the mean of the waytotal-based thresholded fdt_paths. Say thresh = 5%(waytotal), i.e.
So calculate the %waytotal needed in excel then apply to each:
for each subject:
fslstats <subject>/dti/dti.bedpostX/<tractdir>/fdt_paths.nii.gz -l <5%waytotal> -M -V >> meanstreamlines5percent.txt
then import these numbers into Calc and divide by waytotal to get comparable numbers.
#Apply threshold to all pathways to obtain binary masks of the thresholded pathway distribution, to be used to extract values:
then fslmaths <subject>/dti/dti.bedpostX/<tractdir>/fdt_paths.nii.gz -thr <5%waytotal> -bin <subject>/dti/dti.bedpostX/<tractdir>/fdt_paths_thrwaytotal5percent_bin.nii.gz
#Extract Streamlines
fslstats <subject>/dti/dti.bedpostX/<tractdir>/fdt_paths.nii.gz -l <5%waytotal> -M -V >> meanstreamlines5percent.txt
#Now Apply the binarized masks to the FA and radial/axial/mean diffusivity and streamline maps:
#Get Mean FA file
fslstats <subject>/dti/dtifit_FA.nii.gz -k <subject>/dti/dti.bedpostX/<tractdir>/fdt_paths_thrwaytotal5percent_bin.nii.gz -M >> pathway_meanFA.txt
# Calculate radial diffusivity:
foreach s (mody*)
fslmaths $s/dti/dtifit_L2 -add $s/dti/dtifit_L3 -div 2 $s/dti/radialDiff
end
# Get Mean Mean Diffusivity file
foreach subject ( mody028 mody096 mody125 mody175 mody191 mody195 )
fslstats $subject/dti/dtifit_MD.nii.gz -k $subject/dti/dti.bedpostX/probtrackx_lhTPJ2lhsupFrontal/fdt_paths_2percentthr_bin.nii.gz -M >> TPJ2supfrontal-proj0_meanMD.txt
end
# Get Mean Axial Diff file
foreach subject ( mody028 mody096 mody125 mody175 mody191 mody195 )
fslstats $subject/dti/dtifit_L1.nii.gz -k $subject/dti/dti.bedpostX/probtrackx_lhTPJ2lhsupFrontal/fdt_paths_2percentthr_bin.nii.gz -M >> TPJ2supfrontal-proj0_meanAxialDiff.txt
end
In Practice, and consecutively for each subject and each tract of interest, changing the threshold accordingly:
(The key here is to copy-paste the command and use the find and replace function in Word to change what needs to be, then paste the whole page into a terminal) *****************
echo mody028 >> ribbon_streamlines
echo supFrontal >> ribbon_streamlines
fslstats mody028/dti/dti.bedpostX/probtrackx_lhTPJ2lhsupFrontal/fdt_paths.nii.gz -l 1.8 -M -V >> ribbon_streamlines
foreach subject ( mody028 mody096 mody125 mody175 mody191 mody195 )
foreach? fslstats $subject/dti/dtifit_FA.nii.gz -k $subject/dti/dti.bedpostX/probtrackx_lhTPJ2lhsupFrontal/fdt_paths_2percentthr_bin.nii.gz -M >> TPJ2supfrontal-proj0_meanFA.txt
end
and similarly for MD, rD, aD
#TBSS:
mkdir <tbss-dir>
cp <subjects_dtifit_FA><tbss-dir>/<grp>_<alphabetic-name>_FA.nii.gz
cd <tbss-dir>
#check range of data
fslstats <grp>_<alphabetic-name>_FA.nii.gz -R
#should return something between 0 and 1 (roughly) otherwise will need scaling
tbss_1_preproc *
# QA: look at the URL given at the end of the process
# in seychelles: ssh seychelles, nmrenv
cd <tbss-dir>
pbsubmit -c "tbss_2_reg -n" -m cherif@nmr.mgh.harvard.edu -o '-l nodes=1:opteron' #for calculating most representative subject o a fast opteron node
# can come out of seychelles for rest:
tbss_3_postreg -S
cd stats
fslview all_FA mean_FA_skeleton #to QA that the skeleton looks fine: change colormap of skeleton to green (display 2000:8000), turn on movie loop, and see that skeleton falls where it should for all subjects
tbss_4_prestats <thresh> #the threshold for the skeleton image... a common one is 2000... lower it for more tracts, but be careful that they still make sense
#this creates stats/all_FA_skeletonized
# Now for stats, use randomise, for which you need a design matrix design.mat and contrasts design.con
design_ttest2 design <Ngroup1><Ngroup2> # to create these files for a simple 2-group comparison
#otherwise use Glm GUI, but careful that the order of the design entries must match the alphabetical order of the input FAs, which you can check by "ls -l FAi/*FAi.hdr"
randomise -i all_FA_skeletonised -o tbss -m mean_FA_skeleton_mask -d design.mat -t design .con -n <Npermutations, 2000 exploratory, 5000 for real> -c 3 -V
example:
pbsubmit -c "randomise -i all_FA_skeletonised.nii.gz -o tbss_tfce_RT_demean_agecov -m mean_FA_skeleton_mask -d Demeaned_RT_Age.mat -t Demeaned_RT_Age.con -T" -m cherif -o '-l nodes=1:opteron'
pbsubmit -c "randomise -i all_FA_skeletonised.nii.gz -o tbss_tfce_RT_demean_agecov_D -D -m mean_FA_skeleton_mask -d Demeaned_RT_Age.mat -t Demeaned_RT_Age.con -T" -m cherif
pbsubmit -c "randomise -i all_FA_skeletonised.nii.gz -o tbss_tfce_Acc_demean_agecov_D -D -m mean_FA_skeleton_mask -d Demeaned_Acc_age.mat -t Demeaned_Acc_age.con -T" -m cherif
pbsubmit -c "randomise -i all_FA_skeletonised.nii.gz -o tbss_tfce_Acc_demean_agecov -m mean_FA_skeleton_mask -d Demeaned_Acc_age.mat -t Demeaned_Acc_age.con -T" -m cherif
TAKING FOREVER!! SO:
pbsubmit -c "randomise -i all_FA_skeletonised.nii.gz -o tbss_5000_RT_demean_agecov_D -D -m mean_FA_skeleton_mask -d Demeaned_RT_Age.mat -t Demeaned_RT_Age.con -n 5000 -c 3" -m cherif -o '-l nodes=1:opteron'
pbsubmit -c "randomise -i all_FA_skeletonised.nii.gz -o tbss_5000_Acc_demean_agecov_D -D -m mean_FA_skeleton_mask -d Demeaned_Acc_age.mat -t Demeaned_Acc_age.con -n 5000 -c 3" -m cherif -o '-l nodes=1:opteron'
#outputs: tbss_tstat<#> = raw/unthresholded tstats; tbss_maxc_tstat<#> = cluster p-value image (fully corrected for multiple comparison across space - note they are 1-p for convenience of display, so 0.95 is significant)
SO MAXC_TSTAT and MAX_TFCE_TSTAT are the corrected maps, to look at with a threshold of 0.95
#example of display:
fslview MNI152 mean_FA_skeleton -l Green -b 2000,8000 tbss_tstat1 -l Red-Yellow -b 3,6 tbss_tstat2 -l Blue-Lightblue -b 3,6
#Filling the statistic for display purposes:
foreach stat (*5000*maxc*)
foreach? tbss_fill $stat 0.95 mean_FA $stat.fill
foreach? end
To get FA values in areas of sig. diff. from TBSS:
use cluster image to mask all_FA_skeletonised in fslmeants. e.g.
fslmaths your_tfce_max_tstatX_results -thr 0.95 -bin your_results_mask
then
fslmeants -i all_FA_skeletonised -m your_results_mask -o mean_FA_results.txt
TBSS deproject:
If you would like to extract FA values (or anything else in native space for that matter) from the TBSS results:
1) Draw an ROI on the MNI-space result map
2) run tbss_deproject:
tbss_deproject <ROIfile> 1
this will deproject from the skeleton to the all_FA maps
From that you can see where, on the nonlinearly aligned individual FA maps, the blob's voxels came from:
fslview all_FA <ROIfile>_to_all_FA
Alternatively (and more interestingly):
tbss_deproject <ROIfile> 2
this does what the "1" option above does, but ALSO reverts the nonlinear registration, so that now you get one file per subject, <subject>_FA_FAi_<ROIfile>.nii.gz, in the native DTI space.
3) foreach file (*_FA_FAi*)
fslmaths $file -bin $file.bin
4) cd ../FAi
foreach subject (CTRL_mody028 CTRL_mody050 CTRL_mody065 CTRL_mody073 CTRL_mody096 CTRL_mody100 CTRL_mody125 CTRL_mody126 CTRL_mody169 CTRL_mody177 CTRL_mody179 CTRL_mody181 HFA_mody175 HFA_mody187 HFA_mody188 HFA_mody189 HFA_mody192 HFA_mody194 HFA_mody195 HFA_mody197 HFA_mody199 HFA_mody201 HFA_mody202 HFA_mody203)
fslstats ${subject}_FA_FAi -k ../stats/${subject}_FA_FAi_masktest_contrast33_skeleton.nii.gz.bin.nii.gz -M
MORE Simply, on the SKeleton:
foreach contrast (tstat1 tstat2 tstat19 tstat20 tstat21 tstat22 tstat23 tstat24 tstat25 tstat26)
cluster -i 9v12_RT_demean_maxc_${contrast}.nii.gz -t 0.95 -o clusters_${contrast}
fslmaths clusters_${contrast}.nii.gz -thr 1 -uthr 1.9 ${contrast}_cluster1
fslmaths clusters_${contrast}.nii.gz -thr 2 -uthr 2.9 ${contrast}_cluster2
fslmaths clusters_${contrast}.nii.gz -thr 3 -uthr 3.9 ${contrast}_cluster3
fslmaths clusters_${contrast}.nii.gz -thr 4 -uthr 4.9 ${contrast}_cluster4
fslmaths clusters_${contrast}.nii.gz -thr 5 -uthr 5.9 ${contrast}_cluster5
fslmaths clusters_${contrast}.nii.gz -thr 6 -uthr 6.9 ${contrast}_cluster6
fslmaths clusters_${contrast}.nii.gz -thr 7 -uthr 7.9 ${contrast}_cluster7
fslmaths clusters_${contrast}.nii.gz -thr 8 -uthr 8.9 ${contrast}_cluster8
echo ${contrast} >> ${contrast}.txt
foreach cluster (cluster1 cluster2 cluster3 cluster4 cluster5 cluster6 cluster7 cluster8)
echo ${cluster} >> ${contrast}.txt
fslmeants -i all_FA_skeletonised.nii.gz -m ${contrast}_${cluster}.nii.gz >> ${contrast}.txt
end
end