Working with ROIs (Freesurfer)
This page seriously needs to be reorganized. It documents some procedures that were developed before we figured out better ways to do some things, and so it may be more confusing than helpful.
Working From Group-Level Results
Group-level contrast maps, generated by Mri_glmfit are well-suited for generating a set of ROIs that are appropriate for use across the set of participants entered into the group-level analysis. They have the useful property of being spatially consistent across all participants. There are a couple things to note, however. The first is that there must first be a sufficient number of participants to produce statistically significant contrast clusters. The second is that the group-level analysis must be done in fsaverage space. If the data you wish to explore is in native self space, you would need to map the ROIs back to a subject's native surface space. Note that this is not a strong consideration, given that you would have had to have preprocessed and analyzed individuals in fsaverage space in order to produce group-level contrast maps, so I cannot imagine a situation that would absolutely require you to map ROIs to native space.
Working From Existing Labels
The other primary source of ROIs come from anatomical labels provided in the fsaverage/label/ folder. Because the folder in the FreeSurfer install directory is owned by root, you won't be able to add new labels to this folder. Instead, you should make a copy of the fsaverage folder in your project directory:
cp -r $FREESURFER_HOME/subjects/fsaverage $SUBJECTS_DIR
AnchorBar
All the procedures and code lower down on the page have been replaced with a handy Python tool called AnchorBar (aka "Analyzer"), which maintains a sqlite3 database of annotations in template space. The program allows the user to add new annotations to the database and compute set operations on the stored annotations to create new annotation schemes; for example dividing functional clusters along anatomical boundaries.
Adding Annotations to the AnchorBar database
The initial AnchorBar annotation database is called annot.db, and has been populated with standard FreeSurfer annotations.
See also: Working with Subcortical ROIs (Freesurfer)
Old Information
The information below this point is accurate, but is largely obsolete. I couldn't bring myself to just delete it all because you never know when it'll be handy to know some of this stuff.
The .label files created during autorecon3, or by saving an overlay are plaintext lists of vertices. It should be straightforward to find the intersection or union of the vertices in two or more .label files using common Unix shell commands.
.label Syntax
Below is the first few lines of the lh.banksts.label file for a participant. This file was produced during the Lausanne parcellation procedure.
#!ascii label , from subject FS_0043 vox2ras=TkReg 1457 32992 -47.867 -62.891 33.007 0.0000000000 32993 -48.881 -62.907 32.963 0.0000000000 34014 -47.298 -61.696 33.538 0.0000000000 34015 -47.972 -61.259 33.648 0.0000000000 34024 -47.634 -62.264 33.251 0.0000000000 34025 -48.286 -61.775 33.303 0.0000000000 ...etc
The first line is a header line. The second line indicates the number of vertices in the label file. The remaining lines indicate the vertex number, x,y,z coordinates of the vertex, and I have no idea what that last column indicates but it's not terribly important for our purposes.
Working From Individual-Level Results
Running selxavg3-sess will produce some number of statistical contrast maps that can be loaded as an overlay in tksurfer. Significant clusters can be turned into ROIs for that particular individual in a series of steps. This has the advantage of being applicable to a single dataset -- there is no need to wait until sufficient fMRI data has been collected to produce significant group-level contrasts. The downside is that this is process has many steps.
Load Overlay
The first step is to load up the participant's surface data and the desired overlay. Care must be taken in selecting a contrast overlay that is appropriate for the type of analysis that you eventually plan to do. As a rule of thumb, the overlay contrast should be orthogonal (i.e., completely independent) of the contrast that you are ultimately interested in. For example, if your goal is to compare mean signal strength for low- vs. high-familiar trials, it would be inappropriate to use the contrast of high-familiar.vs.low-familiar to generate ROIs from which to calculate the mean signal strength. The reason for this is that this contrast will be selecting voxels that have already been found to demonstrate a significant difference. Naturally, any analysis of familiarity effects on these voxels will be biased! A more appropriate contrast to use would be something like task.vs.rest, which isn't biased towards either high or low-familiarity items, but is instead just showing the voxels that were more highly active during the task trials (which presumably contain equal numbers of high and low familiarity items).
After the overlay is loaded, you may want to play around with the overlay configuration parameters: for example you might choose to truncate your results to only show positive t-values, or modify your significance threshold.
Save Clusters as Labels
Now comes the somewhat tedious part of individually converting each cluster into a .label file. This simple but repetitive process is described here and you can get pretty quick at it:
- If you wish to subdivide a large cluster, first create one or more paths that completely cross (e.g., bisect) the cluster.
- Use the cursor to select a point anywhere in the desired cluster
- Click the Custom Fill button. A new window will popup:
- Choose Up to functional values below threshold
- Optionally choose Up to and including paths if had subdivided the cluster. You will have to repeat these steps on the remainder of the cluster
- The region will be colored in white before reverting to a yellow outlined region
- On the Tools window, go to 'File > Label > Save Selected Label
- A box will open with a path to save the ROI in the subject's label directory by default
- Give it a meaningful name that includes the hemisphere and uses the .label filename extension (e.g. lh.TASK_001.label)
- Click OK
Merge .label Files
If you wish to work with the resulting .label files individually, your work is done. However, if you wish to merge multiple clusters into a single .annot file, you have to first create a color table which is called by mris_label2annot.
Create a CTAB File
Existing instructions for ctab files are a little hard to follow, so I'll do my best. The CTAB files are plaintext files with 6 columns and 1 row for each label. The structure of the columns is as follows:
| Index | Label | R | G | B | A |
You might find it easiest if left and right hemispheres have their own CTAB files. These two files can be identical if you want, but, having two files allows you to have different numbers of labels associated with each hemispheres without encountering problems in subsequent steps. The indices for each hemisphere should be numbered independently. For example, if there are 11 labels in the lh and 8 in the rh, then one ctab file should have indices 1 to 11 and the other should have indices 1 to 8.
The label names should correspond to the ?h.*.label files you created. For example, suppose you created the following 5 label files in the previous step:
lh.TASK_001.label lh.TASK_002.label rh.TASK_003.label rh.TASK_004.label rh.TASK_005.label
The lh ctab file should contain the following labels:
0 unknown 128 128 128 0 1 TASK_001 110 112 100 0 2 TASK_002 128 140 130 0
and the rh ctab file should contain the following labels:
0 unknown 128 128 128 0 1 TASK_003 110 112 100 0 2 TASK_004 190 140 105 0 3 TASK_005 200 150 100 0
Note that as explained above, the ctab indices (column 1) were independent. Sets of RGB values can be generated according to your whims, or you can use an external website to produce a set of n distinguishable colors. For example, I went to http://tools.medialab.sciences-po.fr/iwanthue/ to generate 34 distinct colors. Note that the final column, A, refers to the alpha transparency value, and is always set to 0 as far as I can tell.
Call mris_label2annot
The last step uses the information in the ctab files to generate the .annot files from the appropriate sets of .label files. The syntax is:
mris_label2annot --s ${SUBJECT_ID} \
--h [lh | rh] \
--ctab [lh | rh].${CTAB_FILENAME} \
--a ${PREFIX} \
--ldir ${SUBJECTS_DIR}/${SUBJECT_ID}/label
This will parse the indicated ctab file to identify the labels indicated therein. The identified labels will then be assigned the colors indicated in that ctab file and written to the .annot file specified by the ${PREFIX} associated with the --a argument. It is assumed that all the .label files you created earlier were saved to ${SUBJECTS_DIR}/${SUBJECT_ID}/label. If this is different, you should specify the actual file location with the --ldir argument.
Note that the above command is the syntax used when you have an unknown label and corresponding entry in your CTAB file. If you do not wish to include the unknown label, rather than assign index 0 to the unknown label, start your CTAB index numbering at 0 for the first real label, and use the --no-unknown command-line argument when calling mris_label2annot. I didn't remember this being a problem the first time I tried this, but when I just went through this process again, it seemed that the CTAB index needed to start at 0, but I couldn't have index 0 associated with unknown if there was not an ?h.unknown.label file, even when using the --no-unknown switch. When I edited the ctab files to have index 0, 1, ... , n-1 corresponding with label_001, label_002, ... , label_n, then everything seemed to work fine.
Preparation of .label Files
The group analysis used the fsaverage surface, and since mris_divide_parcellation (which you may also be using, or use later) assumes we're subdividing a particular subject's annotation, it makes sense to copy the source .annot files into the fsaverage label/ folder (giving the appropriate ?h prefix). The group analysis annot file will be found in the contrast directory for a particular analysis. The example below describes working with the .annot file for the TASK contrast in the group analysis, found under $SUBJECTS_DIR/RFX
HEMI=rh
cd $SUBJECTS_DIR/RFX/analysis.${HEMI}/task/glm.wls/osgm
ANNOT_FILE=cache.th30.abs.sign.ocn.annot
cp ${ANNOT_FILE} $SUBJECTS_DIR\fsaverage\label\${HEMI}.${ANNOT_FILE}
Note that fsaverage may not be writable if it is a symbolic link to $FREESURFER_DIR/subjects/fsaverage. If this is the case, you will have to make a copy of fsaverage that you own:
#first unlink the symbolic link to fsaverage cd $SUBJECTS_DIR unlink fsaverage #now make a new empty fsaverage directory. You will own it and so should be able to read/write anything it contains mkdir fsaverage #finally, copy everything in the real fsaverage directory to the directory that you own: cp -r $FREESURFER_DIR/subjects/fsaverage/* fsaverage
Breaking Up .annot into multiple .label Files
The steps described below do not work on .annot files because they are not plaintext and easily manipulated. However, mri_annotation2label will convert a .annot file to a series of .label files:
BASE_LABEL_NAME=TD_W
ANNOT=lh.cache.th30.abs.sign.ocn.annot
HEMI=lh
SURF=orig
mri_annotation2label --annotation ${ANNOT} --hemi ${HEMI} --subject fsaverage --surf ${SURF} --labelbase ${HEMI}.${BASE_LABEL_NAME}
Alternatively, you might want to dump all your .label files into a single output directory. Instead of supplying a labelbase, you specify an output dir (outdir):
mri_annotation2label --subject fsaverage --hemi lh --outdir LDT.lh --annotation lh.ldt_hi.annot
This command created a new subdirectory (LDT.lh), and broke up the ld.ldt_hi.annot into its constituent label files, placing them in the subdirectory.
Note: when working with functional contrast map annotation files, the output files will be called ?h.cluster-nnn.label. If you convert multiple .annot files for the same hemisphere (e.g., from two separate contrasts) with the intention of later merging them, you will run the risk of overwriting earlier .label files unless you take care to rename them as you go.
Set Operations on .label Files
Intersecting .label Files (Intersection)
Given what we know about the contents of a label file, the unix comm command should suffice to find the common entries between two label files, which is the set intersection of vertices appearing in each.
Merging .label Files (Union)
A union operation will find all vertices appearing in either .label file. You will want to make sure that only unique values are retained (i.e., you don't want to list the same vertex multiple times in the output .label file).
The easiest way to accomplish this is FreeSurfer's mri_mergelabels tool:
mri_mergelabels -i label1 -i label2 ... -o outputlabel
or
mri_mergelabels -d <dirname> -o outputlabel
This FreeSurfer tool merges two or more label files. It does this by catting the label files together (after removing the first two lines). It inserts a new header (the first two lines). The number of entries in the new file (ie, the number on the second line), is computed by summing those from the input files. When you pass a directory name with the -d flag, it will merge all labels in a directory. This is handy when you want to merge many labels created from multiple .annot files.
For example, if I have followed the steps above to create a series of ?h.CONTRAST_A-0??.label and ?h.CONTRAST_B-0??.label files, then I can find the union of CONTRAST_A and CONTRAST_B by creating a lh/ and rh/ subdirectory and moving the label files to the appropriate directories. Then:
mri_mergelabels -d lh -o lh.AB_UNION.label mri_mergelabels -d rh -o rh.AB_UNION.label
The above two commands will generate my lh and rh label files in the current working directory (assumed to be ${SUBJECTS_DIR}/${SUBJECT}/label; ${SUBJECT} is likely to be fsaverage if you are making a functional mask from group-level analyses in fsaverage space).
You're Going to Want to Turn These into an .annot File
After you've merged multiple .label files into a single .label file, all the vertices will be assigned to the same label. You will probably want to go through and relabel each of the distinct regions (e.g., by cluster) in FreeSurfer. After you've done that, save (export) each of the relabeled regions as separate .label files (e.g., CLUST_01, CLUST_02, ... etc.) in a working directory.
After you do, you need to create or snag a CLUT file.
Finally, you use mris_label2annot to merge your united labels, along with the CLUT into a single .annot file.
ROI 'VennDiagram' Matlab Function
ROIVennDiagram.m uses intersect and setdiff functions to extract intersecting and unique vertices. Variable names suggest that it's used for T1 and T2 data, but you can use it for any paired data. The label files generated above (using mri_annotation2label) need to be copied into a lh and rh directory within a parent dir. A simple modification would allow these files to remain in your fsaverage or other subject's dir, if you wanted to clutter those up.
function ROIvenndiagram = ROIvenndiagram(hemi)
%% function ROIvenndiagram = ROIvenndiagram('hemi');
%Pull out intersecting and different vertices for FS label files,
%and store in one of 3 text files (T1, T1T2, T2)
%
%Currently takes one string input indicating lh or rh.
%
%You'll need to adjust the file names on lines 42 and 43 to match your
%file names.
%
%Note - The headers in the output text files have quotes around them, as
%the headers begin with matlab unfriendly characters (#). So, you'll need to
%go into those 3 output text files and delete the quotes.
%
%Intended for use in parent dir containing lh and rh dirs, which
%contain label files.
% _____________________________________
% G.J. Smith
% gjsmith4@buffalo.edu
% _____________________________________
parentdir=pwd;
hemidir = sprintf('%s/%s', pwd, hemi);
cd(hemidir);
%% Pre-initialize array for cating labels
T1Agg = cell(1,5);
T2Agg = cell(1,5);
% Begin cycle through label files
for labeln = 1:3000 %arbirary number higher than any label number
%% Apply appropriate number of 0's before number - eg. 001, 010, 100
numcorrection=num2str(labeln);
if labeln < 10
labelnum=strcat('00',numcorrection);
elseif labeln < 100
labelnum=strcat('0',numcorrection);
else
labelnum=strcat(numcorrection);
end
%% T1 and T2 label file names -- Edit these to match your file naming convention
T1_filename = sprintf('%s.T1-%s.label', hemi, labelnum);
T2_filename = sprintf('%s.T2-%s.label', hemi, labelnum);
%% Build one large aggregate file with all vertices from all label files
header = 2; % lines 1 and 2 of label file stored seperately
delimiter = ' ';
if exist([pwd filesep T1_filename], 'file')
%load data
T1file = importdata(T1_filename,delimiter,header);
%Store in one array
T1Agg = vertcat(T1Agg, num2cell(T1file.data));
end
if exist([pwd filesep T2_filename], 'file')
%load data
T2file = importdata(T2_filename,delimiter,header);
%Store in one array
T2Agg = vertcat(T2Agg, num2cell(T2file.data));
end
end
%% Check that label files were found
exist T1file;
if ans == 1
display('Labels found, pulling out intersects and setdiffs');
else
display('No label files found, check your dir or file names');
return;
end
%% Find intersection and differences
T1Agg = cell2mat(T1Agg);
T2Agg = cell2mat(T2Agg);
%compare
%union = union(T1file.data(:,1), T2file.data(:,1));
intersection = intersect(T1Agg, T2Agg, 'rows'); %'rows' works with same # of columns
T1_diffs = setdiff(T1Agg, T2Agg, 'rows');
T2_diffs = setdiff(T2Agg, T1Agg, 'rows');
%% Build new label files
%T1
T1_differences = cell(50000, 5); %hard coded arbitrarily high number of rows for pre-init
T1_differences(1,1) = T1file.textdata(1,1); % header
T1_differences{2,1} = length(T1_diffs); % header line 2, new number of rows
T1_differences(3:length(T1_diffs) + 2,:) = num2cell(T1_diffs); % data
%T2
T2_differences = cell(50000, 5); %hard coded arbitrarily high number of rows for pre-init
T2_differences(1,1) = T2file.textdata(1,1); % header
T2_differences{2,1} = length(T2_diffs); % header line 2, new number of rows
T2_differences(3:length(T2_diffs) + 2,:) = num2cell(T2_diffs); % data
%Match
Match_labelfile = cell(50000, 5); %hard coded arbitrarily high number of rows for pre-init
Match_labelfile(1,1) = T1file.textdata(1,1); % header T1 and T2 for same subject is the same
Match_labelfile{2,1} = length(intersection); % header line 2, new number of rows
Match_labelfile(3:length(intersection) + 2, :) = num2cell(intersection); % data
%% Remove blank rows from oversized pre-alo arrays
T1_differences(all(cellfun(@isempty,T1_differences),2),:) = [];
T2_differences(all(cellfun(@isempty,T2_differences),2),:) = [];
Match_labelfile(all(cellfun(@isempty,Match_labelfile),2),:) = [];
%% Write files
cd(parentdir);
T1_fname = sprintf('%s.T1.label.txt', hemi);
writetable(cell2table(T1_differences), T1_fname,'Delimiter',' ','WriteVariableNames',false);
%T2
T2_fname = sprintf('%s.T2.label.txt', hemi);
writetable(cell2table(T2_differences), T2_fname,'Delimiter',' ','WriteVariableNames',false);
%Match
Match_fname = sprintf('%s.T1T2.label.txt', hemi);
writetable(cell2table(Match_labelfile), Match_fname,'Delimiter',' ','WriteVariableNames',false);
end
Visualize and Cluster (Intersected Label File)
- in progress*
Open tksurfer.
tksurfer fsaverage ?h inflated
Import ?h.aparc.annot files for a reference
Load the intersected label file.
Use cut line to make cluster boundaries (If necessary, otherwise see below).
Erase aparc.annot labels that overlap with a desired label.
Select area you want to cluster and use custom fill to create new label. Up to and including path, up to other labels, and up to unlabeled, filling from last clicked vertex. If you left adjacent aparc.annot labels these can be used as boundaries.
Save each label file using save selected label. -- must premake text files as cant specify file names. -- must be a better way of doing this, or getting export annotation to work.
Pro tip: Colour your label file red, when creating new labels, they will be blue. (see below)
CTAB and Annot
Now you'll need to create a ctab file for the labels to make a .annot file.
If you don't have a ton of labels, this can easily be done by hand by copying and editing the 'FreeSurferColorLUT.txt' file found in /usr/local/freesurfer, and simply renaming the regions(label names) to match your labelling convention. Then save it as a .annot.ctab file, e.g. 'RH.annot.ctab'.
Then you'll want to use the mris_label2annot command found here. There are 2 ways to do this. The first is to pass it each of the label files using the--l switch, and a ctab fie. Example:
mris_label2annot --l rh.label-001.txt --l rh.label-002.txt --l rh.label-003.txt --l rh.label-004.txt --l rh.label-005.txt --l rh.label-006.txt --l rh.label-007.txt --l rh.label-008.txt --l rh.label-009.txt --l rh.label-010.txt --l rh.label-011.txt --l rh.label-012.txt --l rh.label-013.txt --l rh.label-014.txt --l rh.label-015.txt --l rh.label-016.txt --l rh.label-017.txt --l rh.label-018.txt --l rh.label-019.txt --l rh.label-020.txt --l rh.label-021.txt --l rh.label-022.txt --l rh.label-023.txt --l rh.label-024.txt --l rh.label-025.txt --l rh.label-026.txt --ctab RH.annot.ctab --s fsaverage --h rh --annot T1T2intersected
That's unwieldy (and yet, that's how I first did this). The easier way is to throw all the label files into a single directory, and pass the directory and a ctab file.
Now that you have .annot files for your intersecting vertices, you can subparcellate and/or continue with network analyses.
Using R to Work With Colortables
This is kind of a stub topic, but there is an R library when trying to extract a CLUT from an existing .annot file. I introduce the freesurferformats R package. In particular, there is a function to extract the color lookup table from a .annot file:
> library("freesurferformats")
> annotfile="/path/to/subject/label/lh.myannot.annot"
> annot = read.fs.annot(annotfile);
> colortable = colortable.from.annot(annot);
> head(colortable);
If you install this library, you can use this handy R script that I wrote that can be run from the terminal to extract a CLUT from a .annot file.