Sunday, July 1, 2012

Region of Interest Analysis


Before we get down to regions of interest, a few words about the recent heat wave: It's taken a toll. The past few days I've walked out the door and straight into a slow broil, that great yellowish orb pasted in the cloudless sky like a sticker, beating down waves of heat that saps all the energy out of your muscles. You try to get started with a couple miles down the country roads, crenelated spires of heat radiating from the blacktop, a thin rime of salt and dust coating every inch of your skin, and realize the only sound in this inferno is the soles of your shoes slapping the cracked asphalt. Out here, even the dogs have given up barking. Any grass unprotected by the shade of nearby trees has withered and died, entire lawns turned into fields of dry, yellow, lifeless straw. Flensed remains of dogs and cattle and unlucky travelers lie in the street, bones bleached by the sun, eyeless sockets gazing skyward like the expired votaries of some angry sun god.

In short, it's been pretty brutal.

Regions of Interest

Region of Interest (ROI) analysis in neuroimaging refers to selecting a cluster of voxels or brain region a priori (or, also very common, a posteriori) when investigating a region for effects. This can be done either by creating a small search space (typically a sphere with a radius of N voxels), or based on anatomical atlases available through programs like SPM or downloadable from web. ROI analysis has the advantage of mitigating the fiendish multiple comparisons problem, in which a search space of potentially hundreds of thousands of voxels is reduced to a smaller, more tractable area, thus reducing overly stringent multiple comparisons correction thresholds. At first glance this makes sense, given that you may not be interested in a whole brain analysis (i.e., searching for activation in every single voxel in the entire volume); however, it can also be abused to carry out confirmatory analyses after you have already determined where a cluster of activation is.

Simple example: You carry out a whole brain analysis, and find a cluster of fifty voxels extending over the superior frontal sulcus. This is not a large enough cluster extent to pass cluster correction at the whole brain level, but you go ahead anyway and perform an additional ROI analysis focused on the center of the cluster. There are not any real safeguards against this measure, as it is impossible to know what the researcher had in mind when they conducted the test. For instance, what if an investigator happened to simply make a mistake and see the results of a whole brain analysis before doing an ROI analysis? Pretend that he didn't see them? These are questions which may be addressed in a future post about a Bayesian approach to fMRI, but for now, be aware that there exists significant potential for misuse of this technique.

Additional Considerations


Non-Independence
Colloquially known as "double-dipping," non-independence has become an increasingly important issue over the years as ROI analyses have become more common (see Kriegeskorte et al, 2009, 2010).  In order to avoid biasing an ROI toward certain regressors, it is essential that the ROI and the contrast of interest share no common regressors. Consider a hypothetical experiment with three regressors: A, B, and C.  The contrast A-B is used to define an ROI, and the experimenter then decides to test the contrast of A-C within this ROI.  As this ROI is already biased toward voxels that are more active in response to regressor A, this is a biased contrast to conduct. This is not endemic only to fMRI data, but applies to any other statistical comparison where bias is a potential issue.

Correction for Multiple ROI Analyses
Ideally, each ROI analysis should be treated as an independent test, and should be corrected for multiple comparisons. It is better practice to have an a priori ROI that will be used for a single test, instead of exploring several ROIs and then correcting for multiple comparisons afterwards.

ROI Analysis in AFNI

Within AFNI, there exists a useful program called 3dUndump which requires x, y, and z coordinates (in millimeters), radius size of the sphere, and the master dataset where the sphere will be applied. A typical command looks like:

3dUndump -prefix (OutputDataset) -master (MasterDataset) -srad (Radius of Sphere, in mm) -xyz (X, Y, and Z coordinates of sphere)

One thing to keep in mind is the orientation of the master dataset. For example, the standard template that AFNI warps has a positive to negative gradient when going from posterior to anterior; in other words, values in the Y-direction will be negative when moving forward of the anterior commissure. Thus, it is important to note the space and orientation of the coordinates off of which you are basing your ROI, and make sure it matches up with the orientation of the dataset you are applying the ROI to. In short, look at your data after you have generated the ROI to make sure that it looks reasonable.

The following is a short Python wrapper I made for 3dUndump. Those already familiar with using 3dUndump may not find much use in it, but for me, having an interactive prompt is useful:


#!usr/bin/env python

import os
import math
import sys

#Read in Raw user input, assign to variables
print("MakeSpheres.py")
print("Created by Andrew Jahn, Indiana University 03.14.2012")
prefix = raw_input("Please enter the output filename of the sphere: ")
master = raw_input("Please enter name of master dataset (e.g., anat_final.+tlrc): ")
rad = raw_input("Please enter radius of sphere (in mm): ")
xCoord = raw_input("Please enter x coordinate of sphere (MNI): ")
yCoord = raw_input("Please enter y coordinate of sphere (MNI): ")
zCoord = raw_input("Please enter z coordinate of sphere (MNI): ")

#Make string of coordinates (e.g., 0 36 12)
xyzString = xCoord + " " + yCoord + " " + zCoord
printXYZString = 'echo ' + xyzString + ' > sphere_' + rad + 'mm_'+xCoord+'_'+yCoord+'_'+zCoord+'.txt'
os.system(printXYZString) #prints xyzstring to filename given above

#Will need sphere file in this format for makeSpheres function
xyzFile = 'sphere_' + rad + 'mm_'+xCoord+'_'+yCoord+'_'+zCoord+'.txt'

def makeSpheres(prefix, master, rad, xyz ):
cmdString = '3dUndump -prefix '+prefix+ ' -master '+master+' -srad '+rad+' -xyz '+xyz 
os.system(cmdString)
return

makeSpheres(prefix=prefix, master=master, rad=rad, xyz=xyzFile)




Which will generate something like this (Based on a 5mm sphere centered on coordinates 0, 30, 20):


Beta values, time course information, etc., can then be extracted from within this restricted region.

ROI Analysis in SPM (Functional ROIs)

This next example will focus on how to do ROI analysis in SPM through MarsBar, a toolbox available here if you don't already have it installed. In addition to walking through ROI analysis in SPM, this will also serve as a guide to creating functional ROIs. Functional ROIs are based on results from other contrasts or interactions, which ideally should be independent of the test to be investigated within that ROI; else, you run the risk of double-dipping (see the "Non-Independence" section above).

After installation, you should see Marsbar as an option in the SPM toolbox dropdown menu:

1. Extract ROIs

After installing Marsbar, select it from the toolbox dropdown menu.  After Marsbar boots up, click on the menu “ROI Definition”.  Select “Get SPM clusters”.


This will prompt the user to supply an SPM.mat file containing the contrast of interest.  Select the SPM.mat file that you want, and click “Done”.

Select the contrast of interest just as you would when visualizing any other contrast in SPM.  Select the appropriate contrast and the threshold criteria that you want.




When your SPM appears on the template brain, navigate to the cluster you wish to extract, and click “current cluster” from the menu, underneath “p-value”.  Highlight the cluster you want in the SPM Results table.  The highlighted cluster should turn red after you select it.



Navigate back to the SPM results menu, and click on “Write ROIs”.  This option will only be available if you have highlighted a clutter.  Click on “Write one cluster”.  You can enter your description and label for the ROI; leaving these as default is fine.  Select an appropriate name for your ROI, and save it to an appropriate folder.




2. Export ROI

Next, you will need to convert your ROI from a .mat file to an image file (e.g., .nii or .img format).  Select “ROI Definition” from the Marsbar menu and then “Export”.



You will then be presented with an option for where to save the ROI. 



Select the ROI you want to export, and click Done.  Select the directory where you wish to output the image.  You should now see the ROI you exported, with a .nii extension.

3. ROI Analysis

You can now use your saved image for an ROI analysis.  Boot up the Results menu as you would to normally look at a contrast, but when prompted for “ROI Analysis?”, select “Yes”. You can select your ROI from the “Saved File” option; alternatively, you can mask activity based on atlas parcellations by selecting the “Atlas” option.



The constriction of search space will mean fewer multiple comparisons need to be corrected for, and thus increases the statistical power of your contrast.




26 comments:

  1. This comment has been removed by the author.

    ReplyDelete
  2. This comment has been removed by the author.

    ReplyDelete
    Replies
    1. Hi Anannya, can you give some more details? Do these ALFF maps refer to Amplitude of Low Frequency Fluctuations? What exactly are you trying to partition?

      -Andy

      Delete
  3. Sir , yes the ALFF refers to Amplitude of Low Frequency Fluctuations . I want to partition each alff maps into ROI in terms of automated anatomical labeling atlas .

    ReplyDelete
    Replies
    1. How do you want to partition the alff maps into an ROI? What do these maps look like? (I haven't dealt with them before, so I'm not familiar with them.) Do you mean you want to extract some kind of alff parameter from the ROI?

      Delete
  4. Hi Andy, thanks for the blog and videos. It's very helpful.
    How can I do the ROI analysis in SPM12? When I click the result button it does not offer the option of ROI analysis so i don't know how to view the result... I tried matlab help spm_vb_roi_basis but it can not find the ROI files.. any suggestions? thank you

    ReplyDelete
    Replies
    1. Hi Petr,

      The GUI should ask you whether you want to do an ROI analysis, after you specify a contrast. Try doing a simple contrast, like [1 -1], and see whether it then prompts you for doing an ROI analysis.

      Note that the ROI analysis through the GUI is similar to correcting for multiple comparisons within the ROI. If you want to extract data from the ROI, you can either right-click on the glass brain that shows up and select "extract betas," or you can follow the instructions here: http://andysbrainblog.blogspot.com/2014/07/quick-and-efficient-roi-analysis-using.html

      -Andy

      Delete
  5. Hi Andy,
    Thank you very much for all your videos and the great blog!
    In SPM12, if you look at the results of a selected contrast, you have the option to select "Small Volume" which is located under the "whole brain" button. Than you can select "Sphere", "Box", or "Image".
    Do you think it is statistically legal to load mask of an a-priori defined mask (left and right hippocampus, created with the Anatomy toolbox), and preform a small-volume correction like this?

    Best regards,

    Marco

    ReplyDelete
  6. Hey Marco,

    If you have a good reason for using the mask, then go ahead and use it! Yes, it's considered biased to use an SVC once you've already looked at the results, but one would have to prove that you weren't going to use that mask in the first place. And for someone to know what you were thinking would first require access to the mind of God.

    -Andy

    ReplyDelete
  7. Hello Andy! Thank you for the great blog.

    Is there any way to have atlas labels and % overlay from TFCE log(p) images from a VBM analysis?

    best,
    Kostas

    ReplyDelete
    Replies
    1. Hey Kostas,

      Are you referring to depicting your log(p) images on a template brain? You can do that in FSLview, select the Atlas toolbar, and add a region to display on the template. If you overlay both your results and the probabilistic regions, you should be able to see where your results overlap with which atlas labels.

      -Andy

      Delete
  8. Hi Andrew,

    I have some a priori hypotheses about the hippocampus. I am trying to find a way to segment the hippocampus in Marsbar so I can get functionnal values during an encoding task. However, the litterature isn't very clear about HOW to segment the hippocampus. Some researchers segment it manually while others use complex programs to do it, some just arbitrarily put a sphere on the top of the hippocampus etc.

    Do you have any suggestions? I'd like to do it manually using Marsbar, but where to start? What coordinates should I use?

    Thank you in advance, and happy new year!

    Nick

    ReplyDelete
    Replies
    1. Hey Nick,

      To segment a region as small as the hippocampus, you need to make sure that your scans have resolution that's high enough to partition the structure into meaningful units.

      If you have the spatial resolution, then the next thing I would recommend is to look at each hippocampus within the subject's native space, instead of normalizing each subject and extracting data from a standard ROI. This will give you greater precision when identifying what data comes from which subregion.

      That means that the actual segmentation is up to you. When I see study's of this type, usually they mention that they have one or two trained anatomists parcellate the structure on a voxel by voxel basis, and then an ROI analysis is done on those partitions for each subject, and statistics are done across subjects. That's the method I would use.

      -Andy

      Delete
  9. Hey Andrew,

    I have 4 different groups and I do Voxel based analyses (DTI) in SPM: F-contrast and T-contrast: and for testing the interaction, I don't know which combination should I use.

    the 4 different groups: A1, A2, B1, B2 => A

    I used these combination:
    1 1 -1 -1
    1 -1 0 0
    1 0 -1 0
    0 1 0 -1
    0 0 1 -1

    But I'm not sure of these combinations will tell me something about the interaction between these different groups.

    Could you please help me here for!

    Kind regards,

    Zahra

    ReplyDelete
    Replies
    1. Hi Zahra,

      What does the => A mean? Are A1 and A2 subgroups of group A, like different timepoints within subject?

      The typical contrast vector for a 2x2 interaction is [1 -1 -1 1], which will reveal which voxels show a significant interaction. However, to determine what is driving the interaction, I would make four separate contrast vectors, one for each group:

      [1 0 0 0]
      [0 1 0 0]
      [0 0 1 0]
      [0 0 0 1]

      To create contrast images for each group. Then, extract those parameter estimates from an ROI (e.g., a DTI measurement such as FA or MD) and plot them. This will give you an idea whether it is a crossover interaction, or something else.


      Best,

      -Andy

      Delete
  10. Hi Andy, thank you for such insightful videos. Similar problem to above, I'm not prompted on whether id like to complete an ROI analysis in the results box after clicking results. I am using SPM12. Any help would be amazing as I'm trying to complete an analysis for my undergraduate dissertation.

    Many thanks in advance

    ReplyDelete
    Replies
    1. That's odd; after you select a particular contrast, it should prompt you first for whether you want to apply masking, and then whether you want to do an ROI analysis. Have you selected or created a contrast?

      -Andy

      Delete
  11. Hi Andrew,
    I think I have the same questions as above: in SPM 12, after selecting a contrast. It will first prompt whether to apply any mask. But then will be the threshod options: FEW or uncorrect.
    Not sure how to do an ROI for exploration (by saying this I mean ROI before seeing the whole brain activation results? Any suggestions?

    Thank you!
    Pengfei

    ReplyDelete
    Replies
    1. Hi Pengfei,

      Those correction options are within the mask. SPM will do a cluster simulation (or FWE Bonferroni correction) given the restricted search space, and any results you see in the table take into account the correction you applied.

      Best,

      -Andy

      Delete
  12. "... it is good practice to have a priori hypotheses about where activation is to be found in order to avoid these issues of multiple comparisons."
    So if one starts the study with a strong a priori hypothesis of where to locate the ROIs, then a threshold correction based on the number of ROIs does not necessarily apply?
    Thanks
    Geoff

    ReplyDelete
    Replies
    1. Hi Geoff,

      You should still correct for the number of ROI analyses you do for a statistical test; what I meant was, if you have a single ROI that you've planned to use in advance, instead of trying out several ROIs in the hopes that one of them will show a significant effect, and correcting for multiple comparisons after the fact.

      I've changed the wording of that paragraph to be more in line with this thinking.


      Best,

      -Andy

      Delete
  13. Thanks, Andy for the response and for creating such a useful resource.
    Do you find that most papers reporting ROI analyses, don't carry out this correction? I get that impression from looking at a small sample.

    ReplyDelete
    Replies
    1. I haven't paid that much attention to it, but I think it's uncommon. It hadn't crossed my mind until a reviewer asked me to do it.

      -Andy

      Delete
  14. Hi Andy,

    First, thanks a lot for your great blog of functional imaging.

    We have an issue concerning the signal change comparison between different ROIs. Funnily enough you seem to be the only scientist explicitly mentioning it (besides this threat on marsbar's blog (https://sourceforge.net/p/marsbar/mailman/message/31271216/). Do you have any reason (or references anyway) which would speak against using an ANOVA with Region as one of the factors (assuming the other one being for example Treatment–medicine or placebo). The strongest argument against the ANOVA is the physioanatomical differences between different ROIs, which make baseline proportions problematic. Generally, we use a priori selected ROIs as independent tests, and then we Boferroni correct them for the cumulative number of tests.

    For completeness, here is a paper from de Hann (https://www.ncbi.nlm.nih.gov/pubmed/24137123), which used direct comparisons between ROIs, by showing that baselines were not different among the ROIs.

    Thank you very much in advance for any hint on this.

    Emil.

    ReplyDelete
    Replies
    1. Hello Emil,

      Thank you for the reference; I hadn't thought that much about it, aside from the possibility there may be systematic differences between grey matter and white matter physiology. I've heard other researchers talk about certain regions that, due to their proximity to major blood vessels, can show significant differences from other regions, but I haven't read much about this.

      In any case, I believe that using parameter estimates representing the magnitude of the BOLD response is more defensible than percent signal change (which is why in my papers I report parameter and contrast estimates instead of percent signal change). Doing so avoids many of the problems with baseline that are brought up in that email thread. In any case, this is worth keeping in mind, and I'll think about it some more.

      Best,

      -Andy

      Delete
  15. Hendrik HeinbockelJuly 11, 2017 at 7:44 AM

    Hi Andy,

    at first i want to thank you for your videos. They have been a great help to me :).

    I have got a little problem using marsbar for an ROI analysis. I want to extract (high-pass filtered) time-series from my fMRI Data. Though the pipeline is absolutely clear and easy, resulting values to not match what i expect from time-series to be like. Here is an example output from 1 ROI (AAL_Amygdala L) and of 4 Scans:

    443.231
    442.523
    443.562
    440.921

    Plotted, marsbar tells me that those valaues are made of the mean signal intensity, which is indeed what i want.
    But, as i recall time-sereis usually are ranged between about + 3 and -3.
    I would be really thankful for any kind of tip you can provide...

    Best regards,

    Hendrik

    ReplyDelete