#!/bin/tcsh -xef echo "auto-generated by afni_proc.py, Wed Sep 27 11:57:37 2023" echo "(version 7.60, August 21, 2023)" echo "execution started: `date`" # to execute via tcsh: # tcsh -xef proc.FT.simple |& tee output.proc.FT.simple # to execute via bash: # tcsh -xef proc.FT.simple 2>&1 | tee output.proc.FT.simple # =========================== auto block: setup ============================ # script setup # take note of the AFNI version afni -ver # check that the current AFNI version is recent enough afni_history -check_date 14 Nov 2022 if ( $status ) then echo "** this script requires newer AFNI binaries (than 14 Nov 2022)" echo " (consider: @update.afni.binaries -defaults)" exit endif # the user may specify a single subject to run with if ( $#argv > 0 ) then set subj = $argv[1] else set subj = FT.simple endif # assign output directory name set output_dir = $subj.results # verify that the results directory does not yet exist if ( -d $output_dir ) then echo output dir "$subj.results" already exists exit endif # set list of runs set runs = (`count -digits 2 1 3`) # create results and stimuli directories mkdir -p $output_dir mkdir $output_dir/stimuli # copy stim files into stimulus directory cp FT/AV1_vis.txt FT/AV2_aud.txt $output_dir/stimuli # copy anatomy to results dir 3dcopy FT/FT_anat+orig $output_dir/FT_anat # ============================ auto block: tcat ============================ # apply 3dTcat to copy input dsets to results dir, # while removing the first 2 TRs 3dTcat -prefix $output_dir/pb00.$subj.r01.tcat FT/FT_epi_r1+orig'[2..$]' 3dTcat -prefix $output_dir/pb00.$subj.r02.tcat FT/FT_epi_r2+orig'[2..$]' 3dTcat -prefix $output_dir/pb00.$subj.r03.tcat FT/FT_epi_r3+orig'[2..$]' # and make note of repetitions (TRs) per run set tr_counts = ( 150 150 150 ) # ------------------------------------------------------- # enter the results directory (can begin processing data) cd $output_dir # --------------------------------------------------------- # QC: look for columns of high variance find_variance_lines.tcsh -polort 3 -nerode 2 \ -rdir vlines.pb00.tcat \ pb00.$subj.r*.tcat+orig.HEAD |& tee out.vlines.pb00.tcat.txt # ========================== auto block: outcount ========================== # QC: compute outlier fraction for each volume touch out.pre_ss_warn.txt foreach run ( $runs ) 3dToutcount -automask -fraction -polort 3 -legendre \ pb00.$subj.r$run.tcat+orig > outcount.r$run.1D # outliers at TR 0 might suggest pre-steady state TRs if ( `1deval -a outcount.r$run.1D"{0}" -expr "step(a-0.4)"` ) then echo "** TR #0 outliers: possible pre-steady state TRs in run $run" \ >> out.pre_ss_warn.txt endif end # catenate outlier counts into a single time series cat outcount.r*.1D > outcount_rall.1D # ================================= tshift ================================= # time shift data so all slice timing is the same foreach run ( $runs ) 3dTshift -tzero 0 -quintic -prefix pb01.$subj.r$run.tshift \ pb00.$subj.r$run.tcat+orig end # -------------------------------- # extract volreg registration base 3dbucket -prefix vr_base pb01.$subj.r01.tshift+orig"[2]" # ================================= volreg ================================= # align each dset to base volume foreach run ( $runs ) # register each volume to the base image 3dvolreg -verbose -zpad 1 -base vr_base+orig \ -1Dfile dfile.r$run.1D -prefix pb02.$subj.r$run.volreg \ -cubic \ pb01.$subj.r$run.tshift+orig end # make a single file of registration params cat dfile.r*.1D > dfile_rall.1D # compute motion magnitude time series: the Euclidean norm # (sqrt(sum squares)) of the motion parameter derivatives 1d_tool.py -infile dfile_rall.1D -set_nruns 3 \ -derivative -collapse_cols euclidean_norm \ -write motion_${subj}_enorm.1D # create an anat_final dataset, aligned with stats 3dcopy FT_anat+orig anat_final.$subj # ================================== blur ================================== # blur each volume of each run foreach run ( $runs ) 3dmerge -1blur_fwhm 4.0 -doall -prefix pb03.$subj.r$run.blur \ pb02.$subj.r$run.volreg+orig end # ================================== mask ================================== # create 'full_mask' dataset (union mask) foreach run ( $runs ) 3dAutomask -prefix rm.mask_r$run pb03.$subj.r$run.blur+orig end # create union of inputs, output type is byte 3dmask_tool -inputs rm.mask_r*+orig.HEAD -union -prefix full_mask.$subj # ================================= scale ================================== # scale each voxel time series to have a mean of 100 # (be sure no negatives creep in) # (subject to a range of [0,200]) foreach run ( $runs ) 3dTstat -prefix rm.mean_r$run pb03.$subj.r$run.blur+orig 3dcalc -a pb03.$subj.r$run.blur+orig -b rm.mean_r$run+orig \ -expr 'min(200, a/b*100)*step(a)*step(b)' \ -prefix pb04.$subj.r$run.scale end # ================================ regress ================================= # compute de-meaned motion parameters (for use in regression) 1d_tool.py -infile dfile_rall.1D -set_nruns 3 \ -demean -write motion_demean.1D # compute motion parameter derivatives (just to have) 1d_tool.py -infile dfile_rall.1D -set_nruns 3 \ -derivative -demean -write motion_deriv.1D # ------------------------------ # run the regression analysis 3dDeconvolve -input pb04.$subj.r*.scale+orig.HEAD \ -ortvec motion_demean.1D mot_demean \ -polort 3 \ -num_stimts 2 \ -stim_times 1 stimuli/AV1_vis.txt 'BLOCK(20,1)' \ -stim_label 1 Vrel \ -stim_times 2 stimuli/AV2_aud.txt 'BLOCK(20,1)' \ -stim_label 2 Arel \ -gltsym 'SYM: Vrel -Arel' \ -glt_label 1 V-A \ -fout -tout -x1D X.xmat.1D -xjpeg X.jpg \ -fitts fitts.$subj \ -errts errts.${subj} \ -bucket stats.$subj # if 3dDeconvolve fails, terminate the script if ( $status != 0 ) then echo '---------------------------------------' echo '** 3dDeconvolve error, failing...' echo ' (consider the file 3dDeconvolve.err)' exit endif # display any large pairwise correlations from the X-matrix 1d_tool.py -show_cormat_warnings -infile X.xmat.1D |& tee out.cormat_warn.txt # display degrees of freedom info from X-matrix 1d_tool.py -show_df_info -infile X.xmat.1D |& tee out.df_info.txt # create an all_runs dataset to match the fitts, errts, etc. 3dTcat -prefix all_runs.$subj pb04.$subj.r*.scale+orig.HEAD # -------------------------------------------------- # create a temporal signal to noise ratio dataset # signal: if 'scale' block, mean should be 100 # noise : compute standard deviation of errts 3dTstat -mean -prefix rm.signal.all all_runs.$subj+orig 3dTstat -stdev -prefix rm.noise.all errts.${subj}+orig 3dcalc -a rm.signal.all+orig \ -b rm.noise.all+orig \ -expr 'a/b' -prefix TSNR.$subj # --------------------------------------------------- # compute and store GCOR (global correlation average) # (sum of squares of global mean of unit errts) 3dTnorm -norm2 -prefix rm.errts.unit errts.${subj}+orig 3dmaskave -quiet -mask full_mask.$subj+orig rm.errts.unit+orig \ > mean.errts.unit.1D 3dTstat -sos -prefix - mean.errts.unit.1D\' > out.gcor.1D echo "-- GCOR = `cat out.gcor.1D`" # --------------------------------------------------- # compute correlation volume # (per voxel: correlation with masked brain average) 3dmaskave -quiet -mask full_mask.$subj+orig errts.${subj}+orig \ > mean.errts.1D 3dTcorr1D -prefix corr_brain errts.${subj}+orig mean.errts.1D # create ideal files for fixed response stim types 1dcat X.xmat.1D'[12]' > ideal_Vrel.1D 1dcat X.xmat.1D'[13]' > ideal_Arel.1D # -------------------------------------------------- # extract non-baseline regressors from the X-matrix, # then compute their sum 1d_tool.py -infile X.xmat.1D -write_xstim X.stim.xmat.1D 3dTstat -sum -prefix sum_ideal.1D X.stim.xmat.1D # ============================ blur estimation ============================= # compute blur estimates touch blur_est.$subj.1D # start with empty file # create directory for ACF curve files mkdir files_ACF # -- estimate blur for each run in errts -- touch blur.errts.1D # restrict to uncensored TRs, per run foreach run ( $runs ) set trs = `1d_tool.py -infile X.xmat.1D -show_trs_uncensored encoded \ -show_trs_run $run` if ( $trs == "" ) continue 3dFWHMx -detrend -mask full_mask.$subj+orig \ -ACF files_ACF/out.3dFWHMx.ACF.errts.r$run.1D \ errts.${subj}+orig"[$trs]" >> blur.errts.1D end # compute average FWHM blur (from every other row) and append set blurs = ( `3dTstat -mean -prefix - blur.errts.1D'{0..$(2)}'\'` ) echo average errts FWHM blurs: $blurs echo "$blurs # errts FWHM blur estimates" >> blur_est.$subj.1D # compute average ACF blur (from every other row) and append set blurs = ( `3dTstat -mean -prefix - blur.errts.1D'{1..$(2)}'\'` ) echo average errts ACF blurs: $blurs echo "$blurs # errts ACF blur estimates" >> blur_est.$subj.1D # add 3dClustSim results as attributes to any stats dset mkdir files_ClustSim # run Monte Carlo simulations using method 'ACF' set params = ( `grep ACF blur_est.$subj.1D | tail -n 1` ) 3dClustSim -both -mask full_mask.$subj+orig -acf $params[1-3] \ -cmd 3dClustSim.ACF.cmd -prefix files_ClustSim/ClustSim.ACF # run 3drefit to attach 3dClustSim results to stats set cmd = ( `cat 3dClustSim.ACF.cmd` ) $cmd stats.$subj+orig # --------------------------------------------------------- # QC: compute correlations with spherical ~averages @radial_correlate -nfirst 0 -polort 3 -do_clean yes \ -rdir radcor.pb05.regress \ all_runs.$subj+orig.HEAD errts.${subj}+orig.HEAD # ========================= auto block: QC_review ========================== # generate quality control review scripts and HTML report # generate a review script for the unprocessed EPI data gen_epi_review.py -script @epi_review.$subj \ -dsets pb00.$subj.r*.tcat+orig.HEAD # ------------------------------------------------- # generate scripts to review single subject results # (try with defaults, but do not allow bad exit status) # write AP uvars into a simple txt file cat << EOF > out.ap_uvars.txt copy_anat : FT_anat+orig.HEAD mask_dset : full_mask.$subj+orig.HEAD ss_review_dset : out.ss_review.$subj.txt final_epi_dset : vr_base+orig.HEAD vlines_tcat_dir : vlines.pb00.tcat EOF # and convert the txt format to JSON cat out.ap_uvars.txt | afni_python_wrapper.py -eval "data_file_to_json()" \ > out.ap_uvars.json # initialize gen_ss_review_scripts.py with out.ap_uvars.json gen_ss_review_scripts.py -exit0 \ -init_uvars_json out.ap_uvars.json \ -write_uvars_json out.ss_review_uvars.json # ========================== auto block: finalize ========================== # remove temporary files \rm -f rm.* # if the basic subject review script is here, run it # (want this to be the last text output) if ( -e @ss_review_basic ) then ./@ss_review_basic |& tee out.ss_review.$subj.txt # generate html ss review pages # (akin to static images from running @ss_review_driver) apqc_make_tcsh.py -review_style pythonic -subj_dir . \ -uvar_json out.ss_review_uvars.json apqc_make_html.py -qc_dir QC_$subj echo "\nconsider running: \n" echo " afni_open -b $subj.results/QC_$subj/index.html" echo "" endif # return to parent directory (just in case...) cd .. echo "execution finished: `date`" # ========================================================================== # script generated by the command: # # afni_proc.py -subj_id FT.simple -dsets FT/FT_epi_r1+orig.HEAD \ # FT/FT_epi_r2+orig.HEAD FT/FT_epi_r3+orig.HEAD -copy_anat \ # FT/FT_anat+orig -tcat_remove_first_trs 2 -regress_stim_times \ # FT/AV1_vis.txt FT/AV2_aud.txt -regress_stim_labels Vrel Arel \ # -regress_basis 'BLOCK(20,1)' -regress_est_blur_errts -regress_opts_3dD \ # -gltsym 'SYM: Vrel -Arel' -glt_label 1 V-A