3dttest++ -clustsim vs 3dClustSim very different

AFNI Message Board
1

Hi Afni Experts,
We ran 3dttest++ -resid -clustsim on a 2x2x2 MNI grid. The outputs from this stated a minimum cluster of 650 voxels for a p=0.005,alpha=0.05. An error of wrong FD also appeared in the terminal even though the command line processed through to the end.

When running the 3dttest++ residuals with 3dFWHMx -acf and then with 3dClustSim we got 150 voxels for p=0.005, alpha=0.05 with -acf flag selected.

  1. Does 3dttest++ -resid -clustsim automatically use -acf option when calculating smoothness (version downloaded approx 2 weeks ago)?

  2. Is there any reason you could think of as to why the estimates are so far off when the same MNI template mask and residuals are in both scenarios?

Thanks,
Ajay

2

HI Afni Experts,
Regarding question 2 I came across a previous post by Bob which may answer this :
https://afni.nimh.nih.gov/afni/community/board/read.php?1,150994,150996#msg-150996

 I am testing to see if I use 3dTcorr1D with a mask first (mni brainmask x MNI-warped EPI brainmask for overlap), prior to 3dttest++.  This may remove some of the outliers and non-brain regions per subject (along with typically using a mask in 3dttest++).  I will also add the -zskip option since the coverage for each person's MNI warped brain will not be the same (especially cerebellum or part of the cortex).

In addition to my first question in the previous post, Is it recommended to use -resid flag with -clustsim?

Thanks,
Ajay

3

(1) you don’t need or want to use -resid with -Clustsim, since the program will do that for you when you use -Clustsim.

(1) 3dttest++ -Clustsim does NOT use -acf or any smoothness estimates at all. It doesn’t model the spatial structure of the noise directly. It uses a randomization of residuals approach to compute noise-only data, and then runs that stuff through 3dClustsim

(2) The results from 3dttest++ -Clustsim can be very different (and have much larger cluster size thresholds) than 3dFWHMx -acf + 3dClustSim -acf. This is the price for obtaining accurate False Alarm Rate (FAR) control – this is our current approach to dealing with the issues in the Eklund et al. PNAS paper – for more details, see the manuscript http://biorxiv.org/content/early/2016/07/26/065862

(3) I am working on developing a better thresholding tool that still maintains FAR control, but that will be some time off – weeks, at best, considering that even after implementing it, I have to test it, and every test case (as in Eklund et alii and as in our bioRxiv manuscript) takes several days to run – and there are a lot of cases to test. And I have other responsibilities.

4

Hi Bob,
Thanks for the reply. I think after reading another post of yours that the issue may pertain to smoothness estimates of covariates. With the current state of where things are at, do you recommend 3dttest++ -clustsim or using 3dFWHMX -acf on the group level stats (i.e. Deter ended stats Output of 3dMVM if we need it for interactions or residuals from 3dttest++ if that model suffices) along with 3dClustsim for publishing? My concern is that with two groups of 20 and adding age, gender and Moca as covariates that I am losing a lot of power, especially if I’m in the category of FPR of 1-2% of covariates (from your Clustsim post).

Additionally, I have some surface data (i.e. Cortical thickness etc) which I processed through 3dMVM due to the interactions of the model which are more accurate.

  1. Is there any stats tool which directly samples surface data (from my understanding 3dttest++ -clustsim does not work with surfaces just yet)?

  2. If not, is there any surface based smoothness estimation programs through AFNI that use the acf model such as slow_clustsim, or have these not been updated yet from the pure Gaussian model?

Thanks,
Ajay

5

If you are using 3dMVM, I’d say use 3dFWHMx -acf followed by 3dClustSim -acf – there really isn’t another option. However, we don’t have any simulations to back this usage up, since running 1000 3dMVMs on RS data would take approximately forever.

To be clear about the covariates issue in 3dttest++ – the main effect had proper FPR control when the analysis included a covariate, but the covariate effect (slope of the voxel data per subject vs. covariate) did not. That could be fixed – I think – by using the covariate t-statistics in the randomizations that go into 3dClustSim, rather than the main effect t-statistics. If you REALLY need that, I suppose I could provide such an option.

We don’t have any statistics tools worthy of the name for surface-based cluster-size analyses. Considering the complexity of the 3D stuff that exists now and needs to be enhanced significantly (as discussed in the other thread), I don’t see such tools coming into being in AFNI for quite some time – if ever. We aren’t exactly a big group.

6

Hi Bob,
In regards to 3dFWHMx -acf should I use the detrend option and input the stats sub-brick of interest from 3dMVM since no residuals are created? Or is there a better way that I should estimate this?

Thanks for clarifying the covariate issue, I do not think I need this option for the analysis.

As for the surface base I definitely understand it is a huge undertaking for a small group but thanks for giving me the status so I can figure out Plan B!

Best,
Ajay

7

With 3dMVM, there’s a difficulty in that (at present) it does not write out residuals to be analyzed. If you analyze the statistics sub-brick(s), the ACF will be biased by the presence – you hope! – of activations.

One option would be to compute the 3 ACF parameters for each subject’s residuals (errts) time series dataset, then take the mean or median of each of these, and use those 3 values in the 3dClustSim -acf run.

8

Hi Bob,
I just to make sure we are on the same page when you refer to the residuals(errts), for resting state analysis processed through afni_proc.py. Are you referring to the errts.fananticor+orig which went through despike, slice time, volume registration, anaticor, bandpass filter and smoothed (subject space) prior to warping (I had to warp at a separate stage after afni_proc.py due to abnormal brains).

I used 3dFWHMx -detrend -acf ACF_FILE -mask full_mask+orig -input errts.fanaticor+orig. I created an average value for all three parameters which estimates the smoothness as an overall acf smoothness of 4.2 (which is good because I specified 4.0 mm in my afni_proc.py script).

Thanks,
Ajay

9

Hi Ajay,

In the case of resting state analysis with band passing,
it might be better to estimate the smoothness from the
input to the regression, rather than the residuals. The
residuals are what we hope to be BOLD-like, while the
input contains the “noise” that we try to remove.

  • rick
10

Hi Rick,
Since I used 3 PCs for csf and sagittal sinus each, along with anaticor for white matter etc within afni_proc.py, getting all of the “noise” might be a bit difficult as I am not sure I can just add the files together (since some have censored timepoints etc). Would the easiest way to calculate the noise residuals be the following for resting state:

  1. Use 3dTproject to censor pb04 (blurred file) so that it has the same timepoints as errts.fanaticor. This file would represent my signal+noise total file. Then I would substract errts.fanaticor (bold like resting state signal) from pb04.censor.

The result should represent all of the “noise” of csf, sagittal sinus, motion, white matter etc rolled into one file. I would then warp this file to MNI space (same grid as my stats) so that any warping smoothness is added into the estimation. I would then run 3dFWHMx -acf -detrend on this warped noise file to get my estimation of the acf parameters. Finally I would use an average of all subjects to estimate my acf input into 3dClustSim to estimate my clustersize.

Below is the code generated which is run…maybe you can specify which file from there would work if my option above does not seem correct.

  1. On a separate note I noticed my 3dTproject cenmode was zero instead of kill from afni_proc.py. Do you think I should update this or leave as is? I am not sure which is considered the best option currently.

Code Run from AFNI

make ROI PCs and uncensor (zero-pad) : FSvent

3dpc -mask follow_ROI_FSvent+orig -pcsave 3 -prefix roi_pc_01_FSvent
rm.det_pcin_rall+orig"[$ktrs]"
1d_tool.py -censor_fill_parent censor_${subj}_combined_2.1D
-infile roi_pc_01_FSvent_vec.1D -write roi_pc_01_FSvent_noc.1D

make ROI PCs and uncensor (zero-pad) : FSsinus

3dpc -mask follow_ROI_FSsinus+orig -pcsave 3 -prefix roi_pc_02_FSsinus
rm.det_pcin_rall+orig"[$ktrs]"
1d_tool.py -censor_fill_parent censor_${subj}_combined_2.1D
-infile roi_pc_02_FSsinus_vec.1D -write roi_pc_02_FSsinus_noc.1D

------------------------------

run the regression analysis

3dDeconvolve -input pb04.$subj.r*.blur+orig.HEAD
-mask full_mask.$subj+orig
-censor censor_${subj}_combined_2.1D
-ortvec bandpass_rall.1D bandpass
-ortvec roi_pc_01_FSvent_noc.1D ROI.PC.FSvent
-ortvec roi_pc_02_FSsinus_noc.1D ROI.PC.FSsinus
-polort 3 -float
-num_stimts 12
-stim_file 1 motion_demean.1D’[0]’ -stim_base 1 -stim_label 1 roll_01
-stim_file 2 motion_demean.1D’[1]’ -stim_base 2 -stim_label 2 pitch_01
-stim_file 3 motion_demean.1D’[2]’ -stim_base 3 -stim_label 3 yaw_01
-stim_file 4 motion_demean.1D’[3]’ -stim_base 4 -stim_label 4 dS_01
-stim_file 5 motion_demean.1D’[4]’ -stim_base 5 -stim_label 5 dL_01
-stim_file 6 motion_demean.1D’[5]’ -stim_base 6 -stim_label 6 dP_01
-stim_file 7 motion_deriv.1D’[0]’ -stim_base 7 -stim_label 7 roll_02
-stim_file 8 motion_deriv.1D’[1]’ -stim_base 8 -stim_label 8 pitch_02
-stim_file 9 motion_deriv.1D’[2]’ -stim_base 9 -stim_label 9 yaw_02
-stim_file 10 motion_deriv.1D’[3]’ -stim_base 10 -stim_label 10 dS_02
-stim_file 11 motion_deriv.1D’[4]’ -stim_base 11 -stim_label 11 dL_02
-stim_file 12 motion_deriv.1D’[5]’ -stim_base 12 -stim_label 12 dP_02
-fout -tout -x1D X.xmat.1D -xjpeg X.jpg
-x1D_uncensored X.nocensor.xmat.1D
-fitts fitts.$subj
-errts errts.${subj}
-x1D_stop
-bucket stats.$subj

– use 3dTproject to project out regression matrix –

3dTproject -polort 0 -input pb04.$subj.r*.blur+orig.HEAD
-mask full_mask.$subj+orig
-censor censor_${subj}_combined_2.1D -cenmode ZERO
-ort X.nocensor.xmat.1D -prefix errts.${subj}.tproject

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

create an all_runs dataset to match the fitts, errts, etc.

3dTcat -prefix all_runs.$subj pb04.$subj.r*.blur+orig.HEAD

--------------------------------------------------

generate fast ANATICOR result: errts.$subj.fanaticor+orig

--------------------------------------------------

fast ANATICOR: generate local FSWe time series averages

create catenated volreg dataset

3dTcat -prefix rm.all_runs.volreg pb03.$subj.r*.volreg+orig.HEAD

mask white matter before blurring

3dcalc -a rm.all_runs.volreg+orig -b follow_ROI_FSWe+orig
-expr “a*bool(b)” -datum float -prefix rm.all_runs.volreg.mask

generate ANATICOR voxelwise regressors via blur

3dmerge -1blur_fwhm 30 -doall -prefix Local_FSWe_rall
rm.all_runs.volreg.mask+orig

– use 3dTproject to project out regression matrix –

3dTproject -polort 0 -input pb04.$subj.r*.blur+orig.HEAD
-mask full_mask.$subj+orig
-censor censor_${subj}_combined_2.1D -cenmode ZERO
-dsort Local_FSWe_rall+orig
-ort X.nocensor.xmat.1D -prefix errts.$subj.fanaticor

Best,
Ajay

11

Hi Rick,
As a follow up, if you agree with using the pb04.censor-errts.fanaticor option as the best moving foward, would I still need to use 3dFWHMx -detrend option? Or is this no longer needed since there shouldn’t be spatial structure in the noise?

Thanks,
Ajay

12

Hi Ajay,

I am suggesting to use pb04.$subj.r*.blur+orig.HEAD,
which afni_proc.py will demonstrate using the option
-regress_est_blur_epits.

Unfortunately, I still have not put the acf options
in afni_proc.py for those estimations, but you can
add that option.

In any case, yes, use -detrend, and do it one run
at a time.

  • rick
13

Hi Rick,
Thanks for the help. One final question. In my analysis I perform all of my steps in subject space (including creating correlation files from my target seed). I then warp my correlation maps to MNI space. For 3dClustSim would it make sense to calculate my individual ACF parameters for each subject in subject space, determine the average of all subjects, and still provide my MNI mask for the 3dClustSim analysis. My thoughts are that since the pb04 file and analysis is in subject space (and I have full_mask per subject), it may make sense to use these estimates of the noise smoothness. However, I would still supply my MNI mask as the mask since I want to evaluate all voxels within the mask that are in MNI space. The other option is to warp all subjects to MNI and extract the smoothness based on the group mask used in 3dMVM. I just wasn’t sure if it best to include warping in the nosie estimate when the values were calculated in subject space, or if the added smoothness is needed since residuals are “ideally” calculated from the stats residuals (similar to 3dttest++).

Thanks,
Ajay