AFNI version info (afni -ver): Version AFNI_22.3.03 'Lucius Verus'
Hello AFNI message board!
I'm using 3dClustSim to correct the output of a whole brain representational similarity analysis (RSA) I ran in python, but I'm confused about how to generate proper noise estimates for this. Typically statistical significance for these kind of analyses is achieved by permuting behavioural matrices several thousand times to generate null distributions to compare against. This is computationally feasible when run at the resolution of 200-300 ROIs, but not so when run at voxel level as I have here.
For some background - I ran the RSA across two groups that listened to one of two different songs while under the effects of LSD. This outputs a whole brain correlation map for each, which represents the similarity of between-participant a behavioural adjacency matrix, and voxel-wise activity between-participant adjacency matrices.
Having Z-scored these maps, I'm attempting to correct for cluster-size, and then assess which clusters overlap between groups. To do this using 3dClustSim, I understand that I need a noise-estimate that's typically generated by running FWHMx on some distribution of noise whose spatial structure is orthogonal to underlying anatomy.
Currently, I'm generating my noise estimates by running FWHMx on my z-scored RSA outputs with the -2difMAD flag as suggested in the docs, but I don't know whether I should trust this or not, or whether I'm implementing it correctly. What would you do here?
"IF YOUR DATA HAS SMOOTH-ISH SPATIAL STRUCTURE YOU CAN'T GET RID OF:
For example, you only have 1 volume, say from PET imaging. In this case,
the standard estimate of the noise smoothness will be mixed in with the
structure of the background. An approximate way to avoid this problem
is provided with the semi-secret '-2difMAD' option, which uses a combination of
first-neighbor and second-neighbor differences to estimate the smoothness,
rather than just first-neighbor differences, and uses the MAD of the differences
rather than the standard deviation. (If you must know the details, read the
source code in mri_fwhm.c!) [For Jatin Vaidya, March 2010]"
For posterity, I've pasted my current method below.
Thanks in advance for your suggestions!
songs="1 2"
dim="EgoDissolution Simple"
dim="Simple"
wdir="/Users/gcooper/Downloads/FD_040_NoScrub/LSD/rest2/hyperaligned/results"
for d in $dim; do
for song in $songs; do
echo "$song"
echo "$song"
echo "$song"
echo "$song"
echo "$song"
mkdir -p "$wdir"/clustsim_"$song"_"$d"
#z-score maps
MEAN=`3dBrickStat -mean "$wdir"/isrsa_brain_song"$song"_"$d".nii.gz`
STD=`3dBrickStat -stdev "$wdir"/isrsa_brain_song"$song"_"$d".nii.gz`
rm "$wdir"/Z_isrsa_brain_song"$song"_"$d".nii.gz
3dcalc -prefix "$wdir"/Z_isrsa_brain_song"$song"_"$d".nii.gz \
-a "$wdir"/isrsa_brain_song"$song"_"$d".nii.gz \
-expr "(a/a)*((a-$MEAN)/$STD)"
#estimate noise
acf=`3dFWHMx -ACF NULL -2difMAD -input "$wdir"/Z_isrsa_brain_song"$song"_"$d".nii.gz`
acf1=`echo $acf | awk '{print $18}'`
acf2=`echo $acf | awk '{print $19}'`
acf3=`echo $acf | awk '{print $20}'`
echo $acf1 $acf2 $acf3
#run clustsim
3dClustSim -acf $acf1 $acf2 $acf3 \
-nodec -both -iter 10001 \
-pthr 0.05 0.02 0.01 0.005 0.002 0.001 \
-athr 0.05 0.02 0.01 0.005 0.002 0.001 \
-prefix "$wdir"/clustsim_"$song"_"$d"/clustsim_"$song"_"$d"
pthresholds="05 02 01 005 002 001"
cluster_table="$wdir/clustsim_"$song"_"$d"/clustsim_"$song"_"$d".NN2_2sided.1D"
athreshold_col="2"
mkdir "$wdir"/LME_thresh_"$song"_"$d"_a05
#cluster-size correction
for pthreshold in $pthresholds; do
if [ "$pthreshold" = "05" ]
then
cs=`cat $cluster_table | awk -v var="$athreshold_col" 'NR == 9 {print $var}'`
thresh=1.96
echo $cs
fi
if [ "$pthreshold" = "02" ]
then
cs=`cat $cluster_table | awk -v var="$athreshold_col" 'NR == 10 {print $var}'`
thresh=2.33
echo $cs
fi
if [ "$pthreshold" = "01" ]
then
cs=`cat $cluster_table | awk -v var="$athreshold_col" 'NR == 11 {print $var}'`
thresh=2.58
echo $cs
fi
if [ "$pthreshold" = "005" ]
then
cs=`cat $cluster_table | awk -v var="$athreshold_col" 'NR == 12 {print $var}'`
thresh=2.81
echo $cs
fi
if [ "$pthreshold" = "002" ]
then
cs=`cat $cluster_table | awk -v var="$athreshold_col" 'NR == 13 {print $var}'`
thresh=3.09
echo $cs
fi
if [ "$pthreshold" = "001" ]
then
cs=`cat $cluster_table | awk -v var="$athreshold_col" 'NR == 14 {print $var}'`
thresh=3.29
echo $cs
fi
if [ "$pthreshold" = "0005" ]
then
cs=`cat $cluster_table | awk -v var="$athreshold_col" 'NR == 15 {print $var}'`
thresh=3.48
echo $cs
fi
if [ "$pthreshold" = "0002" ]
then
cs=`cat $cluster_table | awk -v var="$athreshold_col" 'NR == 16 {print $var}'`
thresh=3.72
echo $cs
fi
if [ "$pthreshold" = "0001" ]
then
cs=`cat $cluster_table | awk -v var="$athreshold_col" 'NR == 17 {print $var}'`
thresh=3.89
echo $cs
fi
3dmerge \
-dxyz=1 -1clust 1 "$cs" -2clip -100000000 "$thresh" \
-prefix "$wdir"/LME_thresh_"$song"_"$d"_a05/cs"$cs"_t"$thresh".nii.gz \
"$wdir"/Z_isrsa_brain_song"$song"_"$d".nii.gz
done
3dmerge \
-nozero \
-gnzmean \
-prefix "$wdir"/"$song"_cs_all_thresh_all_a05.nii.gz \
"$wdir"/LME_thresh_"$song"_"$d"_a05/cs*_t*.nii.gz
done
#calculate overlapping clusters between songs
rm "$wdir""$d"_overlaps_cthresh_S1S2_a05.nii.gz
3dcalc -prefix "$wdir"/"$d"_overlaps_cthresh_S1S2_a05.nii.gz \
-a "$wdir"/1_cs_all_thresh_all_a05.nii.gz \
-b "$wdir"/2_cs_all_thresh_all_a05.nii.gz \
-expr 'step(a)+(2*step(b))'
done