I have trouble outputting overlap of ROI using whereami -bmask option

AFNI version info (afni -ver):

Precompiled binary linux_ubuntu_16_64: Jan 6 2023 (Version AFNI_23.0.00 'Commodus')

Hello Folks,

So my goal is to identify the ROIs within my mask, and get the overlap of the ROIs and my mask.
My mask was a intersection of two clustered map outputted from two 2dClusterize commands.
Thus I don't have stats on the dataset, just a mask of 0 and 1s.

   3dclust -dxyz=3 -1Dformat -NN2 atlas_inter+tlrc > results_atlas_inter.1D
   whereami -atlas $atlas -coord_file results_atlas_inter.1D'[1,2,3]' -bmask atlas_inter+tlrc -tab > posROIs_atlas_inter.1D

So I get ROI names from within the several atlases i've inputted in the whereami command, but no overlap.

++ Input coordinates orientation set by default rules to RAI
++ In binary mode ...
++    205 voxels in ROI
+++++++ nearby Atlas structures +++++++

Original input data coordinates in TLRC space

Focus point (LPI)                   	Coord.Space 
 -63 mm [L],  -26 mm [P],    9 mm [S]	{TLRC}                                                          
 -64 mm [L],  -27 mm [P],    9 mm [S]	{MNI}                                                           
 -64 mm [U],  -31 mm [P],   14 mm [S]	{MNI_ANAT}                                                      
Atlas          	Within	Label                           	Prob.	Code
Schaefer_Yeo_17n_400	1.0	LH_SomMotB_Aud_2                	 MPM	45 
Schaefer_Yeo_17n_400	2.0	LH_TempPar_3                    	 MPM	197
Schaefer_Yeo_17n_400	6.0	LH_TempPar_2                    	 MPM	196
Schaefer_Yeo_17n_400	6.0	LH_SalVentAttnA_ParOper_1       	 MPM	86 
Schaefer_Yeo_17n_400	7.0	LH_DefaultB_Temp_6              	 MPM	172

What is strange is that I get this overlap when running similar command on one of the map I use in the intersection.

So for illustration, I'll put the 1D files outputted by 3dclust or 3dclusterize for comparison

First, the outputted file that gives me overlap using wherami :

#
#  Cluster report 
#[ Dataset prefix      = TtestTvsR_BOLD2 ]
#[ Threshold vol       = [1] 'taskVsRest_B_Zscr' ]
#[ Supplement dat vol  = [0] 'taskVsRest_B_mean' ]
#[ Option summary      = 1sided,RIGHT_TAIL,2.45726,clust_nvox,154,NN2 ]
#[ Threshold value(s)  = right-tail thr=2.457260 ]
#[ Aux. stat. info.    =  ]
#[ Nvoxel threshold    = 154;  Volume threshold = 4158.000 ]
#[ Single voxel volume = 27.000 (microliters) ]
#[ Neighbor type, NN   = 2 ]
#[ Voxel datum type    = float ]
#[ Voxel dimensions    = 3.000 mm X 3.000 mm X 3.000 mm ]
#[ Coordinates Order   = RAI ]
#[ Mask                = /home/mococo/Documents/Simon_Boylan2/hick_entropy_analysis/results/meanGMmask+tlrc ]
#[ Mean and SEM based on signed voxel intensities ]
#
#Volume  CM RL  CM AP  CM IS  minRL  maxRL  minAP  maxAP  minIS  maxIS    Mean     SEM    Max Int  MI RL  MI AP  MI IS
#------  -----  -----  -----  -----  -----  -----  -----  -----  -----  -------  -------  -------  -----  -----  -----
   1270   54.9   23.6   19.5   31.5   70.5  -28.5   55.5   -7.5   61.5   7.7758   0.1063   23.053   64.5   22.5   13.5 
    837  -60.1   24.2    7.4  -73.5  -37.5   -4.5   49.5  -10.5   28.5   7.7446   0.1074   17.868  -67.5   13.5    7.5 
    447  -49.2   29.6   53.3  -67.5  -28.5    7.5   61.5   31.5   76.5   7.4763   0.1211   15.161  -43.5   16.5   64.5 
    314  -47.6  -20.9    4.5  -61.5  -31.5  -52.5   -4.5  -10.5   22.5   5.6334   0.1005   11.438  -58.5  -13.5    7.5 
    278   -0.3   54.8  -19.6  -31.5   28.5   37.5   70.5  -37.5    1.5   2.9871   0.0567   8.1115   -1.5   49.5   -1.5 
    201    1.6   -2.4   59.8  -13.5   13.5  -13.5    7.5   46.5   76.5   7.5276   0.2343   17.647    1.5   -1.5   58.5 
#------  -----  -----  -----  -----  -----  -----  -----  -----  -----  -------  -------  -------  -----  -----  -----
#  3347   -3.9   20.7   21.0                                             7.1144   0.0591                             

Now, the broken one that does not :

#
#Cluster report for file /home/simon/Documents/hick_entropy/results/task_versus_rest_clusterization_analysis/Model_results/atlas_inter+tlrc.HEAD 
#[Connectivity radius = 1.44 mm  Volume threshold = 54.00 ]
#[Single voxel volume = 27.0 (microliters) ]
#[Voxel datum type    = byte ]
#[Voxel dimensions    = 3.000 mm X 3.000 mm X 3.000 mm ]
#[Coordinates Order   = RAI ]
#[Fake voxel dimen    = 1.000 mm X 1.000 mm X 1.000 mm ]
#Mean and SEM based on Absolute Value of voxel intensities: 
#
#Volume  CM RL  CM AP  CM IS  minRL  maxRL  minAP  maxAP  minIS  maxIS    Mean     SEM    Max Int  MI RL  MI AP  MI IS
#------  -----  -----  -----  -----  -----  -----  -----  -----  -----  -------  -------  -------  -----  -----  -----
    158   63.3   25.5    9.2   46.5   70.5    4.5   52.5   -4.5   22.5        1        0        1   67.5   10.5   -4.5 
     47    6.2   -0.5   64.4    1.5   10.5  -10.5    7.5   55.5   76.5        1        0        1   10.5    1.5   55.5 
#------  -----  -----  -----  -----  -----  -----  -----  -----  -----  -------  -------  -------  -----  -----  -----
#   205   50.2   19.6   21.9                                                  1        0                             

Can anyone understand why doesn't wherami outputs overlap ?
I didn't find any precisions on -bmask usage in help files.

Thanks

Howdy-

I'm a little confused about using both -bmask .. and -coord_file .. within a single whereami call. I would probably do those separately. In the first case, you are asking about ROI overlap with the reference -atlas .. from a volumetric dataset. With the second opt, you are asking about where the particular point (an "x,y,z" location) falls, within/near the reference atlas. I guess it is fine to do those as a single call---the program doesn't complain---but those are separate things. (Oddly, the program complained about using -dxyz=3 with 3dclust when I tried to run it; note that 3dClusterize is a much more modern, clear and flexible program, and I would highly recommend using that, instead.)

That being said, I tried what you were running on my system, and I did get a list of percentage overlaps of the volumetric cluster map. Looking at the point-location output info, that point does not fall within any particular ROI---the "within" column value is >0.0 in all cases. Are you sure that your atlas_inter+tlrc ROI falls within part of the Schaefer-Yeo atlas of interest? One way to check that could be:

3dcalc -a atlas_inter+tlrc -b DSET_SY -expr 'step(a)*step(b)' -prefix CHECK_OLAP.nii.gz
3dinfo -max CHECK_OLAP.nii.gz

... and see if the max is nonzero. Or:

3dABoverlap \
    -no_automask \
    atlas_inter+tlrc  \
   DSET_SY

... and see what the overlap report text shows.

In terms of reporting overlaps, I think it is probably much nicer to use ROI overlaps than just peak-voxel or center-of-mass values. The single-voxel values can be very misleading and leave out a lot of information. This is discussed more here:

-> in particular in this section.

We do have the program gen_cluster_table that will make tables like the one in Table 2 there, which might be helpful. Note that if the mask region is muuuch bigger than the ROIs of interest, asking for 10% overlap as a minimal amount is probably not appropriate, and that could be heavily lowered.

--pt

Woah,
Thank you so much for all these info.
I had to write a little code myself to do half what 3dABoverlap does.

Results show that 93.68 % of my atlas_inter is within Schaefer_Yeo atlas.
They do not have the same size and I had to resample the atlas to check :

#A=./atlas_inter+tlrc.BRIK  B=/home/simon/Documents/hick_entropy/results/resampleSchaefer_17N_400.nii.gz
#A           #B           #(A uni B)   #(A int B)   #(A \ B)     #(B \ A)     %(A \ B)    %(B \ A)    Rx(B/A)    Ry(B/A)    Rz(B/A)
205          42846        42861        190          15           42656         7.3171     99.5566     2.3557     8.0868     1.5762

Thank you also for the tips on reporting results, this will be very helpful, you don't even imagine.

I must say, I really like AFNI for its flexibility, but it is really hard to have a broad view of the possibilities that already existing. I discover new programs every week.
We would need a lot more tutorials I guess, or to know where to find them.

Still, when I try to run

whereami -atlas Schaefer_Yeo_17n_400 -bmask atlas_inter+tlrc -tab > out.1D

I get this out.1D :

++ Input coordinates orientation set by default rules to RAI
++ In binary mode ...
++    205 voxels in ROI

I will give a look to the links you provided.

Reading your publication, I am pleased that this is actually what I have been doing as I found the voxelwise approach uninformative not easy to read and prone to some biases indeed.

This is actually the reason why I would love to have a quick and easy way to use whereami, with one atlas and my mask.

Thank you very much for your help (again)

1 Like

Sure, I understand about the tutorials and greater ease of "how to" do things.

We are trying to write more demos for things, esp. for papers and for new/larger tools.

  • The AFNI Codex contains AFNI commands and in some cases full sets of commented scripts for processing.
  • We have some tutorials here.

Probably the best place to go for "what program might I use to do this thing" is the Classified Program List. The idea is that you can type keywords to search within the page, or by section and browse the descriptions.

That being said, we can always do more. It is a tricky thing.

--pt

Hello,

I didn't mean to complain ^^.
AFNI is awesome !!!! I would love to participate in the future !

Right now, I'm just lacking resources. So thank you very much for the links.

Any ideas why I don't have any overlap calculated with this code ?

If this is the same dataset that does overlap with the SY atlas, then I will defer to @dglen for insight.

--pt

I'm not sure where the problem might be. Check these things:

AFNI_SUPP_ATLAS_DIR is set to the path location of these atlas files. That can be set in your .afnirc file in your home directory, on the command line (@AfniEnv -set), export/setenv in a shell script or set for each program use on the command line.

For the bmask and omask options, make sure your input data matches the resolution of the atlas. You may want to resample these datasets to match the template and then compute the overlap

3dresample -master ~/abin/MNI152_2009b_template.nii.gz -input mydset.nii.gz -prefix mydset_rsMNI.nii.gz -rmode NN

Still, you should get some report or some error as here:

whereami -DAFNI_SUPP_ATLAS_DIR=`pwd` -bmask testROI12+tlrc. -atlas Schaefer_Yeo_17N_400 -tab
++ Input coordinates orientation set by default rules to RAI
++ In binary mode ...
++    5536 voxels in ROI
++    5536 voxels in atlas-resampled mask
Intersection of ROI (all non-zero values) with atlas Schaefer_Yeo_17N_400 (sb0):
   4.5  % overlap with RH_VisCent_ExStr_9, code 210
   4.3  % overlap with RH_ContA_PFCl_1, code 330
   3.8  % overlap with RH_VisCent_ExStr_11, code 212
   3.4  % overlap with RH_SalVentAttnB_PFCl_1, code 306
   2.9  % overlap with RH_SalVentAttnB_PFCl_3, code 308
   2.8  % overlap with RH_SomMotB_Cent_2, code 257
...

Hello,

I finally found the solution (I think).
I figured out that inter_masks where constructed from two different maps and by default, 3dcalc puts the data in TLRC instead of the native MNI.

This is why whereami didn't find anything, it seems like it assumes that spaces where different.

This is probably unwanted behavior, or at least some warning would be nice.
Some programs run a warning and continues anyway.

The solution if anyone is looking for it is to refit your bmask in the same space as your atlas
refit -space MNI atlas_inter+tlrc.HEAD

Thank you for your help

Hello again,

I managed to use gen_cluster_table function, very handy.
But the help exemple found here is wrong :

Examples ~1~

1) Basic usage to create a table:
   gen_cluster_table                                  \
        -input_map   Clust_mask+tlrc.HEAD             \
        -input_atlas MNI_Glasser_HCP_v1.0.nii.gz      \
        -prefix      table_cluster_olap_glasser.txt

2) Basic usage to create a table, using a lower overlap fraction cut-off:
   gen_cluster_table                                  \
        -input_map     Clust_mask+tlrc.HEAD           \
        -input_atlas   MNI_Glasser_HCP_v1.0.nii.gz    \
        -min_olap_perc 5                              \
        -prefix        table_cluster_olap _glasser.txt

I think it's '-input_clust' now not input_map.

Thank you for your help

Hi-

Well, actually, effectively all datasets marked as being in any reference space get a "+tlrc" AFNI view before .HEAD and .BRIK. For example, even if the dataset were in Haskins pediatric or NIMH Macaque Template (NMT) space, it would still be "+tlrc" in the filename. 3dcalc does not change space or grid.

What was your reference template space from processing? And what is the result of this on your original mask:

3dinfo -space -av_space -prefix DSET

?

-pt

Thanks for pointing on the help file mistake; I will correct that now.

--pt

Hello,

I checked and indeed, half of my datasets are in MNI space but specified as tlrc in the header.
I didn't check how or why, I know they are in the correct space.

Would it be possible to add the space warning in whereami ?
I think some afni functions has it by default.

Thank you for your help.

To be clear, you are saying that if you run:

3dinfo -space -prefix DSET1 DSET2 DSET3 ....

... there is a mix of spaces reported?

If that is the case for data that were meant to go through a single pipeline, that is a pretty major thing to check on. Do you have the analysis scripts still, hopefully?

If the analysis was done in AFNI, you can check the processing history of a file with either:

# general, see all info including history at bottom of terminal text
3dinfo DSET

# specific, view just history
3dinfo -history DSET

The history accumulates with each AFNI program run, so you can track back. Note that running a non-AFNI program will purge the history (it's stored in the AFNI header, and other software don't propagate that, which is reasonable). But again, this is a major thing to be sure of.

I think I see what you mean for whereami: if a dset has space=MNI and someone is asking for a non-MNI-space atlas, then pop out a warning? That is probably a good idea, yes. Trickiness comes in with using atlases that are made by other software that don't have space name info in the same way, but at least in cases where one expects them to match, one would be able to verify.

--pt

Thank you for your concern,

I had the space warning at the beginning of the pipeline, but I checked that it was a false alert, just didn't clear all the datasets that were in TLRC. I think it happened because I had to use tedana, and I think the output was in tlrc for some reason, or something similar. Nevertheless, I checked and all my datasets are in MNI.

Thank you again