error with mean/abs in RBA (Region-Based Analysis Program through Bayesian Multilevel Modeling)

AFNI version info (afni -ver):

Precompiled binary linux_openmp_64: Nov 12 2020 (Version AFNI_20.3.02 'Vespasian')
and/or
Precompiled binary linux_ubuntu_16_64: Jan 6 2023 (Version AFNI_23.0.00 'Commodus')

R version : "R version 4.3.1 (2023-06-16)" Beagle Scouts

Hello folks,

I hope you can help me.
I'm running a group analysis on a tapping task with 26 subjects.
I wanted to compare results from 3dttest++ and RBA (bayesian analysis) but I cannot run RBA without hitting on an error :

Warning message:
In mean.default(lop$dataTable$Y) :
  argument is not numeric or logical: returning NA
Error in abs(r) : non-numeric argument to mathematical function
Calls: fisher -> ifelse
Execution halted

What would it be due to ?

Could you share the RBA script and a few lines of the data table?

Gang

Hello Gang Chen,
Thanks for your help and your work ;)

Here is the script :
RBA -prefix myResult -chains 4 -iterations 1000 -model 1 -EOI 'Intercept'
-r2z -dataTable myData.txt \

And MyData.txt is like such :
Subj ROI CMRO2
AlPu LprimMotCrtx 0.100674
AlPu RprimMotCrtx 0.013181
AlPu LprimVisCrtx 161.413
AlPu RprimVisCrtx 0.064862
AlPu LprimSensCrtx 0.141907
AlPu RprimSensCrtx 0.129364
AlPu LprimAudioCrtx -0.0165933
AlPu RprimAudioCrtx -0.0400383
.... .... ....
SoSe RprimVisCrtx 0.263449
SoSe LprimSensCrtx 807.848
SoSe RprimSensCrtx 2.30982
SoSe LprimAudioCrtx -6.08336
SoSe RprimAudioCrtx -0.22615

Add "-Y CMRO2" to your script and check if RBA functions correctly.

Gang

Thanks Gang for your advice.
It got one step further, but I got a different error this time :

Error in ifelse(abs(r) < 0.995, 0.5 * (log(1 + r) - log(1 - r)), stop("Are you sure that you have correlation values so close to 1 or -1?")) : 
  Are you sure that you have correlation values so close to 1 or -1?
Calls: fisher -> ifelse
In addition: Warning messages:
1: In log(1 + r) : NaNs produced
2: In log(1 - r) : NaNs produced
Execution halted

What does this error suggests ?

The -r2z option is intended to signify that the response variable CMRO2 is in correlation form. However, upon reviewing the provided lines, it appears that the last column does not strictly contain correlation values, as evidenced by numbers as large as 807.848 and -6.08336. Verify the accuracy of these values. If they do not represent correlation values, remove the -r2z option.

Gang

Well indeed, These are not correlation values, but CMRO2 (oxygen).
I took this line from the simplest example of the -help page.
I hadn't notice it was done for corelation values.

Meybe it would be interresting to clarify this point in the help page ?

Thank you for your advices.

I tried without the -r2z option and it is working fine
Thank you

How does the RBA result compare to the Student's t-test?

In the next AFNI release the option -r2z will be removed from the example script.

Gang

Hello Gang,
I selected some ROIs to run RBA, but I can't run the 3dClusterize function after 3dttest, it seems the _Zsc sub-brick of the 3dttest output is not recognize as stat :

++ Looks like clustering is *not* being performed on a stat value; won't worry about sidedness, then.
** FATAL ERROR: Thresholding volume does *not* appear to be a stat so you can't use the 'p=..' internal conversion to stat

By the way, is there a usage of RBA that would be doing a clustering like 3dclustsim and 3dclusterize might ?
Meaning making RBA pick the ROIs that are significant rather than picking them myself ?
That would be a low info PRIOR.
And are the "GLM" results provided by RBA done in a 3dttest-like program ?
If this is the case, something strange is happening : here is what I get at the end of the RBA results text file :

===== Summary of region effects under GLM for Intercept (for reference only): no adjustment for multiplicity =====
                     mean       SD  2-sided-p      2.5%         5%        50%      95%     97.5%
LprimAudioCrtx  1.4396300 2.319816 0.54097764 -28.03643 -2.5362368  1.4396300 16.08637  6.238535
LprimMotCrtx   14.3842557 3.279475 0.00021504 -27.28543  8.7636563 14.3842557 35.09005 21.168367
LprimSensCrtx  21.5073800 3.269197 0.00000103 -20.03171 15.9043957 21.5073800 42.14828 28.270230
LprimVisCrtx    5.7954706 2.643320 0.03872387 -27.79110  1.2651591  5.7954706 22.48474 11.263595
RprimAudioCrtx  0.9365258 2.458326 0.70672852 -30.29946 -3.2767287  0.9365258 16.45778  6.021960
RprimMotCrtx   20.2356883 3.933777 0.00003265 -29.74769 13.4936992 20.2356883 45.07258 28.373327
RprimSensCrtx  22.7944383 3.855129 0.00000500 -26.18962 16.1872427 22.7944383 47.13476 30.769380
RprimVisCrtx    4.3596727 2.919706 0.14897996 -32.73871 -0.6443289  4.3596727 22.79397 10.399546

The thing that strikes me is that the values for the 97.5 percentil are inferior to those of the 95 percentil. That shouldn't be right ?

If I ignore the 2.5 and 97.5 percentils of the CI, results are similar, except for one ROI, where it is ambiguous and I need the 2.5 % lower CI bound to be sure.

I can't run the 3dClusterize function after 3dttest, it seems the _Zsc sub-brick of the 3dttest output is not recognize as stat

What does 3dinfo reveal about the 3dttest++ output?

3dinfo -verb 3dttest++output

Stringent clustering may not always be a valid methodology for result reporting. Firstly, clusters often lean towards being overly conservative due to the questionable assumption inherent in mass univariate modeling. Secondly, their boundaries can be arbitrary, lacking direct neurological relevance. Thirdly, a cluster can span multiple anatomical regions, making it problematic to rely solely on a peak voxel to represent the entire cluster, resulting in significant information loss.

As an alternative approach, one may set a threshold (e.g., a voxel-level p-value of 0.01 or 0.02) and then determine an appropriate cluster threshold based on voxel resolution and anatomical structures (e.g., 20). Subsequently, adopting a "highlight, but don't hide" strategy when presenting results can be more aligned with open science principles and enhance reproducibility, as argued in papers such as this and this. This video may also help.

is there a usage of RBA that would be doing a clustering like 3dclustsim and 3dclusterize might? Meaning making RBA pick the ROIs that are significant rather than picking them myself ?

RBA-based modeling, as explained here, is designed to address multiplicity through an integrative approach; therefore, no additional step is required. Given that the statistical evidence is on a continuous spectrum, you may consider it in conjunction with your domain knowledge.

are the "GLM" results provided by RBA done in a 3dttest-like program ?

The main RBA results in the output should be your focus. The GLM output comes from a separate model at each ROI. In other words, it is provided as a byproduct for reference only.

Gang

1 Like

Hello Gang, and thank you very much for your help.

regarding the GLM, I noticed an error in the percentiles, the 2.5 % and 97.5 % values are aberrant. I.e :

                    ...............   mean       SD     2-sided-p    2.5%       5%        50%        95%     97.5%
LprimAudioCrtx  1.4396300 2.319816 0.54097764 -28.03643 -2.5362368  1.4396300 16.08637  6.238535

See how the 2.5% value is above the 5% and the 97.5% is under the 95% ? There is something wrong in the display maybe.

Regarding my non-stat sub-brick problem, here is the output I have on the output of 3dttest++ :

Number of values stored at each pixel = 2
  -- At sub-brick #0 'Tapping-rest_mean' datum type is float:     -14.5056 to       108.082
  -- At sub-brick #1 'Tapping-rest_Zscr' datum type is float:     -3.66587 to       6.54178
     statcode = fizt

I can use Pval on the sub-brick [1], so it seems it is a p-val with all the data needed.
I am running this function :

3dPval -prefix pval_tapping-rest_BOLD tapping-rest_BOLD+tlrc
3dClusterize -inset 'pval_tapping-rest_BOLD+tlrc' -NN 2 -clust_nvox 129 -pref_map BOLDClusters_cond1 -mask $dataFolder/ArDC/mask_group+tlrc -ithr 1 -idat 0 -2sided  p=0.005

Any clues ?

I will take a look at the documentation you sent me, hopefully it will give me a better way to analyse our data.

P.S : After fideling with AFNI's GUI, I found the t-score that matches the p=0.005 and managed to do the clusterization.
the function p2dsetstat helps getting the t-score value needed.

Did you use the option -toz in 3dttest++? I wonder 3dClusterize would not complain if you rerun your 3dttest++ script without -toz.

The GLM output table does appear to be unusual. I may need to ask you to share the data and RBA script for further diagnosis.

Gang