About brainloc

The brainloc package serves to describe a location on a brain surface mesh that is in MNI305 space, like the surfaces of the fsaverage standard brain templates used by FreeSurfer. Typical use cases of this package in the field of computational neuroimaging include:

Statistical results are often presented as significant clusters of activation on a brain surface mesh (sets of connected vertices). Typical use cases for working with clusters are:

This vignette describes how to perform these tasks with the brainloc package.

Preparation: Getting the data required to run the examples

You will need a directory with the brain templates you want to use, typically these are the fsaverage subjects that come with the FreeSurfer neuroimaging software suite. There are 2 ways to get them:

Using the template subjects from a local FreeSurfer installation

When FreeSurfer is installed and configured properly on your system, there is nothing you need to do. In that case, the environment variable SUBJECTS_DIR points at the subjects directory under the installation directory (FREESURFER_HOME), and it will be used to find the templates by the get_subjects_dir function. You can also set SUBJECTS_DIR to your custom template directory, of course, and that will be used instead.

If you cannot or do not want to set environment variables, you can set the option brainloc.subjects_dir in R to any directory that contains the expected directory structure like this: options("brainloc.subjects_dir"="~/my_fs_home/subjects/"). This would mean that ~/my_fs_home/subjects/ exists on your system and holds the directories of the template subjects. We assume that you have the fsaverage subject in there. Setting this option takes precedence over any environment variables.

Have the brainloc package download them for you from the internet

This requires that you have an internet connection and accept the FreeSurfer software license. See the code in the first example of this vignette below, which downloads the templates. The templates will be downloaded only if they cannot be found, and downloaded ones will be found on the next call to the function, of course.

Usage Examples Part I: Location information for a vertex or coordinate

Creating a brainparc object

Let's say we have a coordinate in MNI305 surface space and want to find the closest vertex. We first load brainloc and create a brainparc object, which internally is a named list holding the surfaces (the left and right hemisphere meshes) and the surface parcellations (according to some brain atlas) for these meshes. The parcellations assign to each vertex of a mesh a single region.

One can construct a brainparc manually from arbitrary loaded meshes and parcellations using the brainparc function, but here we use the function brainparc_fs, which is more convenient if you want to read it from files in a SUBJECTS_DIR of recon-all output data. We create a brainparc for the fsaverage template:

library("brainloc");
# Searches for FreeSurfer SUBJECTS_DIR on machine and downloads if needed.
sjd = get_subjects_dir(allow_download = TRUE, accept_freesurfer_license = TRUE); 
bp = brainparc_fs(sjd, "fsaverage", surface = "white", atlas = "aparc");

By default, many functions in the brainloc package print information to STDOUT, which is handy for creating command line clients from the library. We do not want that here, so we disable it globally here:

options("brainloc.silent"=TRUE);

Note: One can also decide for each function that prints something, they all have a silent parameter.

Finding the vertex closest to a query coordinate

Let's say we have some query coordinates in surface space and want to find the closest vertex in the brainparc for each of the coordinates. Here is how we can do that:

query_coords = matrix(seq(9)+0.1, ncol = 3, byrow = TRUE);
query_coords

ccv = coord_closest_vertex(query_coords, get_surface(bp));
ccv

This assumes that the query coordinates are in the surfaces space and uses Euclidean distance. It works for all subjects and meshes.

Finding the MNI152 and Talairach coordinates for a point or vertex in MNI305 space

This assumes that you query vertex (or coordinate) is in MNI305 space. This is true for all FreeSurfer standard template subjects, like fsaverage, fsaverage6, etc. Our brainparc is for fsaverage, so we can use the coordinates of its vertices to query.

query_vertices = c(10L, 145029L);
query_hemis = c("lh", "rh"); # vertex 10 on left hemi and 145029L on right hemi.
query_coords = get_surface_coords(bp, query_vertices, query_hemis);
query_coords

coord_info = coord_MNI305_info(query_coords, surface = get_surface(bp));
coord_info

Note: If you have the regfusionr package installed, have a look at the method parameter of the coord_MNI305_info function: you can use the more accurate regfusionr method to transform from MNI305 to MNI152 then. Open the docs with ?coord_MNI305_info to learn more.

Finding the brain atlas region a vertex belongs to, and the distance to all other brain atlas regions from the vertex

We can easily find the regions a vertex is assigned to in all atlases of a brain parcellation:

query_vertices = c(10L, 145029L);
query_hemis = c("lh", "rh"); # vertex 10 on left hemi and 145029L on right hemi.
regions = vertex_regions(bp, query_vertices, query_hemis);

regions

One can also compute the distance of a vertex to all atlas regions. Several definitions for the distance of a vertex to a region are available, here we go with the defaults:

query_vertices = c(10L, 145029L);
query_hemis = c("lh", "rh"); # vertex 10 on left hemi and 145029L on right hemi.
region_dists = vertex_closest_regions(bp, query_vertices, query_hemis, num_regions_to_report=5L);

region_dists

One can change the distance definition with the parameters distance and linkage of the vertex_closest_regions function. Open the docs with ?vertex_closest_regions to learn more.

Usage Examples Part II: Location information for a cluster

A cluster is a set of vertices on a surface, and the statistical values that are assigned to the cluster's vertices. Clusters are used in vertex-wise analyses to define activated brain regions (where the exact meaning of activated depends on the context). A clusters has a freeform label or name, but typically only consecutive numbers are used.

The data structure used to store clusters in brainloc is called clusterinfo. The following two concepts are relevant to understanding how a clusterinfo instance is typically constructed:

A clusterinfo object can be constructed for a subject in two ways:

In the following example, we construct a clusterinfo instance from a cluster overlay and the statistical map:

lh_tmap_file = system.file("extdata", "lh.tmap.mgh", package = "brainloc", mustWork = TRUE);
rh_tmap_file = system.file("extdata", "rh.tmap.mgh", package = "brainloc", mustWork = TRUE);
lh_overlay_file = system.file("extdata", "lh.cluster.overlayID.mgh", package = "brainloc", mustWork = TRUE);
rh_overlay_file = system.file("extdata", "rh.cluster.overlayID.mgh", package = "brainloc", mustWork = TRUE);
clinfo = clusterinfo(lh_overlay_file, rh_overlay_file, lh_tmap_file, rh_tmap_file, template_subject = "fsaverage", subjects_dir = sjd);
is.clusterinfo(clinfo);

One can extract the clusters from the `clusterinfo`` instance like this:

clusters = get_clusters(clinfo);

This returns a named list, where the keys are cluster names and the values are integer vectors containing the vertices that belong the the respective cluster. See ?get_clusters for more options, like returning only the clusters for a specific hemisphere.

The clusterinfo object can now be used to compute various information on the clusters. We will start by computing information on the location of the clusters:

cl_details = cluster_location_details(clinfo);
cl_details;

Note that the clusterinfo object contains a brain parcellation (a brainparc object, see at the top). We can also compute the overlap of the cluster with its atlas regions:

overlap = cluster_region_overlap(clinfo);
overlap;

Sometimes one wants to find all peaks of a cluster. A peak exists at vertices which have the most extreme value in their neighborhood. For a brain mesh, the neighborhood of a vertex is defined as all vertices which are reachable by one hop along the mesh edges.

peaks = cluster_peaks(clinfo);
peaks

This concludes the package vignette, thanks for reading.



dfsp-spirit/brainloc documentation built on Jan. 28, 2022, 12:25 p.m.