Skip to contents

What this vignette covers

This guide focuses on loading surface parcellations as actual meshes (neurosurf::LabeledNeuroSurface)—not ggseg flats. We cover:

  • Discovering available surface templates in TemplateFlow
  • Schaefer on fsaverage6 (packaged geometry)
  • Schaefer on fsaverage5/fsaverage (via TemplateFlow, when available)
  • Glasser MMP1.0 (via TemplateFlow, when available)
  • Loading raw surface geometry from TemplateFlow
  • Downloading surfaces to a local folder
  • What files are fetched, what the data slot contains, and how labels map.
  • Minimal patterns to attach your own per-vertex values.

Discovering available surface templates

Before loading surfaces, you can query what’s available in TemplateFlow:

# Find surface-related template spaces
tflow_spaces(pattern = "fs")
#> [1] "fsLR"      "fsaverage"

# List available fsLR surfaces
tflow_files("fsLR", query_args = list(hemi = "L"))

# Filter by specific density and surface type
tflow_files("fsLR", query_args = list(
  hemi = "L",
  density = "32k",
  suffix = "midthickness"
))

# List fsaverage surfaces at a specific resolution (e.g., fsaverage6)
tflow_files("fsaverage", query_args = list(
  hemi = "L",
  resolution = "06"  # fsaverage6
))

Common query parameters for surfaces: - hemi: “L” or “R” - density: “32k”, “164k” (for fsLR) - resolution: “06” (fsaverage6), “05” (fsaverage5) - suffix: “pial”, “white”, “inflated”, “midthickness”, “sphere”

Schaefer surface atlases

# 200 parcels, 7 networks, fsaverage6 inflated (packaged geometry, no TF needed)
atl <- schaefer_surf(
  parcels  = 200,
  networks = 7,
  space    = "fsaverage6",
  surf     = "inflated"
)

Returned structure:

  • atl$lh_atlas, atl$rh_atlas: LabeledNeuroSurface objects.
  • slot(..., "data"): integer parcel IDs per vertex.
  • atl$labels / atl$orig_labels: name lookups; atl$network: network per ID.

Using TemplateFlow spaces (when available)

Note: TemplateFlow surface support requires the fsaverage template to be available in your TemplateFlow installation. As of this writing, surface templates may have limited availability. Use fsaverage6 (above) for reliable access to Schaefer atlases.

# fsaverage5 white surface; uses TemplateFlow geometry (if available)
atl_tf <- schaefer_surf(
  parcels  = 400,
  networks = 17,
  space    = "fsaverage5",
  surf     = "white"
)

If TemplateFlow surfaces are available, this uses the same labels/data layout as above.

Attaching your own values

vals <- rnorm(length(atl$ids))        # one value per parcel
mapped <- map_atlas(atl, vals)        # returns LabeledNeuroSurface per hemi with data filled

map_atlas() writes your parcel values into the data slot of each hemisphere surface.

Glasser MMP1.0 surface atlas (when available)

Note: Like Schaefer TemplateFlow surfaces, Glasser surface support requires the fsaverage template in TemplateFlow. Check availability before using.

# Requires fsaverage template in TemplateFlow
glas <- glasser_surf(space = "fsaverage", surf = "pial")

slot(glas$lh_atlas, "data")[1:5]  # parcel IDs
glas$labels[1:5]                  # names for those IDs

If available, geometry comes from TemplateFlow (fsaverage pial/white/inflated/midthickness); labels come from the HCP-MMP1.0 fsaverage annotations published on Figshare as the “HCP-MMP1_0 projected on fsaverage” dataset (ID 3498446; DOI: 10.6084/m9.figshare.3498446). The data vector length equals the vertex count of the mesh; each entry is the parcel ID at that vertex.

What’s in the data slot?

  • Schaefer / Glasser surface atlases: integer parcel IDs per vertex.
  • If you construct a NeuroSurface yourself (e.g., from load_surface_template()), you supply the per-vertex numeric values.

Saving / reusing

saveRDS(atl, "schaefer200_7_fs6_inflated.rds")
# reload later
atl2 <- readRDS("schaefer200_7_fs6_inflated.rds")

Sanity check: parcel surface renders

We snapshot a labeled Schaefer surface to verify geometry+labels load. Runs only when TemplateFlow is available (for non-packaged spaces).

dir.create("figures", showWarnings = FALSE)
png_path <- file.path("figures", "schaefer200_fs6_L_inflated.png")

atl <- schaefer_surf(
  parcels  = 200,
  networks = 7,
  space    = "fsaverage6",  # packaged geometry (no TF), keeps it fast
  surf     = "inflated"
)

geom_l <- neurosurf::geometry(atl$lh_atlas)
neurosurf::snapshot_surface(geom_l, file = png_path)
png_path

If the image renders, the parcellated surface is displayable.

Loading raw surface geometry from TemplateFlow

To load just the surface geometry (without parcellation labels), use get_surface_template() or load_surface_template():

# Get path to a surface file
fslr_path <- get_surface_template(
  template_id = "fsLR",
  surface_type = "midthickness",
  hemi = "L",
  density = "32k"
)
print(fslr_path)

# Load as a neurosurf SurfaceGeometry object
fslr_geom <- load_surface_template(
  template_id = "fsLR",
  surface_type = "inflated",
  hemi = "L",
  density = "32k"
)

# Load both hemispheres at once
both_hemi <- load_surface_template(
  template_id = "fsLR",
  surface_type = "pial",
  hemi = "both",
  density = "32k"
)
# Returns list with $L and $R

Downloading surfaces to a local folder

To copy TemplateFlow surfaces to a local directory:

# Get paths to surfaces (automatically downloaded to cache)
files <- tflow_files("fsLR", query_args = list(hemi = "L", density = "32k"))

# Copy to local folder
dest_folder <- "~/my_surfaces/fsLR"
dir.create(dest_folder, recursive = TRUE, showWarnings = FALSE)
file.copy(files, dest_folder)

# Or use a helper function for bulk downloads
download_surfaces <- function(space, query_args, dest_folder) {
  dir.create(dest_folder, recursive = TRUE, showWarnings = FALSE)
  files <- tflow_files(space, query_args = query_args)
  if (length(files) > 0) {
    dest_paths <- file.path(dest_folder, basename(files))
    file.copy(files, dest_paths, overwrite = TRUE)
    message("Copied ", length(files), " files to ", dest_folder)
  }
  invisible(dest_paths)
}

# Download all fsLR 32k surfaces
download_surfaces("fsLR", list(density = "32k"), "~/my_surfaces/fsLR_32k")

Common workflows

  • Vertex-wise stats: build your own NeuroSurface from load_surface_template() and write stats into data.
  • Parcel summaries: use volumetric data with reduce_atlas() or surface atlases with map_atlas() to project parcel statistics back onto the mesh for visualization.

Dense vertex-wise overlays

plot_brain() supports dense (per-vertex) overlays in addition to parcel-level colouring. You can pass a list with lh and rh numeric vectors — one value per vertex — and the result is a continuous colour map draped over the cortical surface.

atl <- schaefer_surf(200, 7, space = "fsaverage6", surf = "inflated")

make_spatial_overlay <- function(atlas_hemi) {
  verts <- neurosurf::vertices(neurosurf::geometry(atlas_hemi))
  y <- verts[, 2]
  z <- verts[, 3]
  vals <- sin(y / 25) * cos(z / 30)
  as.numeric(vals)
}

ov <- list(
  lh = make_spatial_overlay(atl$lh_atlas),
  rh = make_spatial_overlay(atl$rh_atlas)
)

plot_brain(
  atl,
  overlay         = ov,
  overlay_alpha   = 0.7,
  overlay_palette = "vik",
  interactive     = FALSE
)

Automatic NeuroVol projection

When overlay is a NeuroVol (e.g. a whole-brain statistical map), plot_brain() automatically projects it onto the surface via neurosurf::vol_to_surf() — no manual preprocessing required. Here we build a synthetic volume with two positive “activation” blobs over frontal cortex and one negative blob over occipital cortex.

spacing <- c(4, 4, 4)
origin  <- c(-80, -120, -60)
dims    <- c(40L, 50L, 40L)

sp  <- neuroim2::NeuroSpace(dims, spacing = spacing, origin = origin)
arr <- array(0, dim = dims)

for (i in seq_len(dims[1])) {
  for (j in seq_len(dims[2])) {
    for (k in seq_len(dims[3])) {
      coord <- origin + (c(i, j, k) - 1) * spacing
      d1 <- sum((coord - c(-20, 30, 40))^2)
      d2 <- sum((coord - c( 20, 30, 40))^2)
      d3 <- sum((coord - c(  0,-60, 20))^2)
      arr[i, j, k] <- 3 * exp(-d1 / 800) +
                       3 * exp(-d2 / 800) -
                       2.5 * exp(-d3 / 600)
    }
  }
}

stat_vol <- neuroim2::NeuroVol(arr, sp)

plot_brain(
  atl,
  overlay           = stat_vol,
  overlay_threshold = 0.5,
  overlay_palette   = "vik",
  overlay_alpha     = 0.7,
  interactive       = FALSE
)

In practice you would load a real statistical map instead:

stat_vol <- neuroim2::read_vol("my_tstat.nii.gz")
plot_brain(atl, overlay = stat_vol, overlay_threshold = 2,
           overlay_palette = "vik", interactive = FALSE)

Two optional parameters control the volume-to-surface projection:

  • overlay_fun: interpolation function — "avg" (default), "nn" (nearest-neighbour), or "mode" (most frequent label).
  • overlay_sampling: sampling strategy — "midpoint" (default, halfway between white and pial), "normal_line" (multiple samples along the surface normal), or "thickness" (cortical-thickness–aware sampling).
plot_brain(
  atl,
  overlay          = stat_vol,
  overlay_fun      = "nn",
  overlay_sampling = "normal_line",
  overlay_threshold = 2,
  interactive       = FALSE
)

Notes on dependencies

  • fsaverage6 geometry is bundled; other spaces require TemplateFlow.
  • Chunk eval = FALSE by default; set to TRUE locally with network and TemplateFlow.