The example script below will run FSL's FDR calculation on a 3rd level cope image. The output will be a p-map (activation map of P-values) masked by the FDR threshold calculated (based on the degrees-of-freedom in the 3rd level analysis), as well as text output file that will list the various thresholds that were calculated. This script is based on information taken from the FDR User Guide.
This script will loop through all copes sending each one to a different computing node for speed.
#!/bin/bash FDRRATE=0.05 for i in 1 2 3 4 5 6 7 8 9 do qsub -v EXPERIMENT=Dummy.01 qsub_fsl_fdr.sh $i $FDRRATE done
Example usage:
[user@head ~]$ qsub -v EXPERIMENT=Dummy.01 qsub_fsl_fdr.sh 1 0.05
Copy and edit script below:
#!/bin/sh # This is a BIAC template script for jobs on the cluster # You have to provide the Experiment on command line # when you submit the job the cluster. # # > qsub -v EXPERIMENT=Dummy.01 script.sh args # # There are 2 USER sections # 1. USER DIRECTIVE: If you want mail notifications when # your job is completed or fails you need to set the # correct email address. # # 2. USER SCRIPT: Add the user script in this section. # Within this section you can access your experiment # folder using $EXPERIMENT. All paths are relative to this variable # eg: $EXPERIMENT/Data $EXPERIMENT/Analysis # By default all terminal output is routed to the " Analysis " # folder under the Experiment directory i.e. $EXPERIMENT/Analysis # To change this path, set the OUTDIR variable in this section # to another location under your experiment folder # eg: OUTDIR=$EXPERIMENT/Analysis/GridOut # By default on successful completion the job will return 0 # If you need to set another return code, set the RETURNCODE # variable in this section. To avoid conflict with system return # codes, set a RETURNCODE higher than 100. # eg: RETURNCODE=110 # Arguments to the USER SCRIPT are accessible in the usual fashion # eg: $1 $2 $3 # The remaining sections are setup related and don't require # modifications for most scripts. They are critical for access # to your data # --- BEGIN GLOBAL DIRECTIVE -- #$ -S /bin/sh #$ -o $HOME/$JOB_NAME.$JOB_ID.out #$ -e $HOME/$JOB_NAME.$JOB_ID.out #$ -m ea # -- END GLOBAL DIRECTIVE -- # -- BEGIN PRE-USER -- #Name of experiment whose data you want to access EXPERIMENT=${EXPERIMENT:?"Experiment not provided"} EXPERIMENT=`findexp $EXPERIMENT` EXPERIMENT=${EXPERIMENT:?"Returned NULL Experiment"} if [ $EXPERIMENT = "ERROR" ] then exit 32 else #Timestamp echo "----JOB [$JOB_NAME.$JOB_ID] START [`date`] on HOST [$HOSTNAME]----" # -- END PRE-USER -- # ********************************************************** # -- BEGIN USER DIRECTIVE -- # Send notifications to the following address #$ -M yourname@somewhere.edu # -- END USER DIRECTIVE -- # -- BEGIN USER SCRIPT -- # User script goes here COPE=$1 FDRRATE=$2 STATSDIR=$EXPERIMENT/Analysis/FSL/3rdLevel/GroupTrialtypeA.gfeat/cope1.feat/stats OUTDIR=$STATSDIR # Convert T-values to P-values using the cope and variance cope images cd $STATSDIR ttologp -logpout logp${COPE} varcope${COPE} cope${COPE} `cat dof` # Caclulate the p-value image fslmaths logp${COPE} -exp p${COPE} # Calculate the FDR from the p-value image fdr -i p${COPE} -m ../mask -q $FDRRATE > fdr${COPE} # Get the FDR threshold echo "FDR (rate) is" $FDRRATE FDRTHRESH=`awk 'NR==2' fdr${COPE}` echo "FDR threshold is" $FDRTHRESH FDRTHRESH=`awk -v a=$FDRTHRESH -v b="1" 'BEGIN{print (b - a)}'` echo "1 minus FDR threshold is" $FDRTHRESH # Calculate a p-value image using the FDR threshold fslmaths p${COPE} -mul -1 -add 1 -thr $FDRTHRESH -mas ../mask thresh_fdr_1_minus_p${COPE} # Clean up unneccesary files rm -f logp${COPE}.nii.gz p${COPE}.nii.gz fdr${COPE} # -- END USER SCRIPT -- # # ********************************************************** # -- BEGIN POST-USER -- echo "----JOB [$JOB_NAME.$JOB_ID] STOP [`date`]----" OUTDIR=${OUTDIR:-$EXPERIMENT/Analysis} mv $HOME/$JOB_NAME.$JOB_ID.out $OUTDIR/$JOB_NAME.cope${COPE}.$JOB_ID.out RETURNCODE=${RETURNCODE:-0} exit $RETURNCODE fi # -- END POST USER--