freesurfer.qmd

```{r lockfile, include=FALSE}
renv::use(
  base64enc = "base64enc@0.1-3",
  bit = "bit@4.5.0",
  bit64 = "bit64@4.5.2",
  bslib = "bslib@0.9.0",
  cachem = "cachem@1.1.0",
  cli = "cli@3.6.4",
  clipr = "clipr@0.8.0",
  colorspace = "colorspace@2.1-1",
  cpp11 = "cpp11@0.5.0",
  crayon = "crayon@1.5.3",
  digest = "digest@0.6.37",
  dplyr = "dplyr@1.1.4",
  evaluate = "evaluate@1.0.3",
  fansi = "fansi@1.0.6",
  farver = "farver@2.1.2",
  fastmap = "fastmap@1.2.0",
  fontawesome = "fontawesome@0.5.3",
  forcats = "forcats@1.0.0",
  fs = "fs@1.6.5",
  generics = "generics@0.1.3",
  GGally = "GGally@2.2.1",
  ggplot2 = "ggplot2@3.5.1",
  ggstats = "ggstats@0.8.0",
  glue = "glue@1.8.0",
  gtable = "gtable@0.3.6",
  highr = "highr@0.11",
  hms = "hms@1.1.3",
  htmltools = "htmltools@0.5.8.1",
  isoband = "isoband@0.2.7",
  jquerylib = "jquerylib@0.1.4",
  jsonlite = "jsonlite@1.9.1",
  knitr = "knitr@1.50",
  labeling = "labeling@0.4.3",
  lattice = "lattice@0.22-6",
  lifecycle = "lifecycle@1.0.4",
  magrittr = "magrittr@2.0.3",
  MASS = "MASS@7.3-60.2",
  Matrix = "Matrix@1.7-0",
  memoise = "memoise@2.0.1",
  mgcv = "mgcv@1.9-1",
  mime = "mime@0.13",
  munsell = "munsell@0.5.1",
  nlme = "nlme@3.1-164",
  patchwork = "patchwork@1.3.0",
  pillar = "pillar@1.11.0",
  pkgconfig = "pkgconfig@2.0.3",
  plyr = "plyr@1.8.9",
  prettyunits = "prettyunits@1.2.0",
  progress = "progress@1.2.3",
  purrr = "purrr@1.0.2",
  R6 = "R6@2.6.1",
  rappdirs = "rappdirs@0.3.3",
  RColorBrewer = "RColorBrewer@1.1-3",
  Rcpp = "Rcpp@1.1.0",
  readr = "readr@2.1.5",
  renv = "renv@1.1.5",
  rlang = "rlang@1.1.5",
  rmarkdown = "rmarkdown@2.29",
  sass = "sass@0.4.9",
  scales = "scales@1.3.0",
  stringi = "stringi@1.8.4",
  stringr = "stringr@1.5.1",
  tibble = "tibble@3.2.1",
  tidyr = "tidyr@1.3.1",
  tidyselect = "tidyselect@1.2.1",
  tinytex = "tinytex@0.56",
  tzdb = "tzdb@0.4.0",
  utf8 = "utf8@1.2.4",
  vctrs = "vctrs@0.6.5",
  viridisLite = "viridisLite@0.4.2",
  vroom = "vroom@1.6.5",
  withr = "withr@3.0.2",
  xfun = "xfun@0.51",
  yaml = "yaml@2.3.10"
)
```

# FreeSurfer Measures {#sec-freesurfer}

![Cortical Mesh (source: https://freesurfer.net)](https://freesurfer.net/landing/images/brain.png)

One of the core scans in the A2CPS neuroimaging protocol is a T1w scan, which is used to capture structural information. These scans are processed using a software suite called [FreeSurfer](https://freesurfer.net/). FreeSurfer takes the 3D image, segments it into pre-defined anatomical regions (like the hippocampus or amygdala), and also creates a triangular "mesh" of cerebral cortex from which it derives measurements of structural morphometry, including cortical thickness, surface area, and curvature. These measurements are related to a wide range of diseases, including chronic pain, and are often used as a source of features for predictive modeling. A2CPS has compiled the measurements into tabular files, providing a rich set of neuroanatomical variables ready for statistical analysis to explore brain-behavior relationships or group differences. This kit provides an overview of those tabular files. 

## Starting Project

### Locate data

In the release, data are stored underneath the `mris/derivatives` folder:

```{bash}
#| eval: false
/corral-secure/projects/A2CPS/products/consortium-data/pre-surgery/mris/derivatives/freesurfer
```

This folder contains all of the FreeSurfer outputs for every participant (e.g., this is the directory that could correspond to the FreeSurfer environment variable `$SUBJECTS_DIR`), which subject folders of the form `sub-[recordid]_ses-[protocolid]`.

Most users will not need the raw FreeSurfer outputs and can instead rely on tables in which the morphological measures have been aggregated. There are three relevant tables

- aparc.tsv 
  - measurements derived from parcellations of the cortical surface (e.g., cortical thickness)
  - produced by [`aparcstats2table`](https://surfer.nmr.mgh.harvard.edu/fswiki/aparcstats2table)
- aseg.tsv 
  - measurements derived from segmentations of the 3d structural image (e.g., average signal intensity)
  - produced by [`asegstats2table`](https://surfer.nmr.mgh.harvard.edu/fswiki/asegstats2table)
- headers.tsv 
  - information about the brain as a whole (e.g., estimated Total Intracranial Volume)
  - derived from the headers of the `aparcstats2table` and `asegstats2table` outputs

For each of these, there are data dictionaries in `json` files (e.g., `aparc.json` explains the columns of `aparc.tsv`). 

```{bash}
#| eval: false
$ ls /corral-secure/projects/A2CPS/products/consortium-data/pre-surgery/mris/derivatives/freesurfer/*{json,tsv}
/corral-secure/projects/A2CPS/products/consortium-data/pre-surgery/mris/derivatives/freesurfer/aparc.json  /corral-secure/projects/A2CPS/products/consortium-data/pre-surgery/mris/derivatives/freesurfer/gm_morph.tsv
/corral-secure/projects/A2CPS/products/consortium-data/pre-surgery/mris/derivatives/freesurfer/aparc.tsv   /corral-secure/projects/A2CPS/products/consortium-data/pre-surgery/mris/derivatives/freesurfer/headers.json
/corral-secure/projects/A2CPS/products/consortium-data/pre-surgery/mris/derivatives/freesurfer/aseg.json   /corral-secure/projects/A2CPS/products/consortium-data/pre-surgery/mris/derivatives/freesurfer/headers.tsv
/corral-secure/projects/A2CPS/products/consortium-data/pre-surgery/mris/derivatives/freesurfer/aseg.tsv
```

Copies of the data dictionaries are available online for preview: [aparc.json](https://github.com/a2cps/snapshot/blob/v2.1.0/src/snapshot/data/aparc.json), [aseg.json](https://github.com/a2cps/snapshot/blob/v2.1.0/src/snapshot/data/aseg.json), [headers.json](https://github.com/a2cps/snapshot/blob/v2.1.0/src/snapshot/data/headers.json).

### Extract data


In this kit, we will explore the cortical parcellations. 

```{r}
library(readr)
library(dplyr)
library(tidyr)
library(ggplot2)
```


As described above, the parcellation information is stored in `aparc.tsv`. 

```{r load}
aparc <- read_tsv("data/aparc.tsv", show_col_types = FALSE)
head(aparc)
```

That table contains nine regional measurements (GrayVol, SurfArea) for the various structures that compose each of six parcellations.

A typical analysis will only rely on one parcellation. For details of the different parcellations, see the [FreeSurfer documentation](https://surfer.nmr.mgh.harvard.edu/fswiki/CorticalParcellation). Two good starting choices are either `aparc` or `aparc.a2009s`. The `aparc` atlas is a coarser parcellation, whereas `aparc.a2009s` is finer. Here, we'll filter for `aparc.a2009s`.

```{r a2009s}
a2009s <- aparc |>
  filter(parc == "aparc.a2009s")
head(a2009s)
```


### Example Analysis: Relationship between Cortical Surface Area and Gray Matter Volume

One of the most common features in predictive models is cortical thickness, but each of the morphological measurements may carry unique information. For example, cortical surface area may be more predictive of widespread chronic pain as compared to cortical thickness, which may be more closely related to chronic headaches [@bhatt_mapping_2024]. Let's compare the regional volume with the regional surface area.

:::  {#fig-pairs}

```{r pairs}
a2009s |>
  ggplot(aes(x = SurfArea, y = GrayVol)) +
  geom_point(alpha = 0.1)

```

Cortical Gray Matter Volume (mm^3) and Surface Area (mm^2). Each point represents a single pair of measurements for a single region within a participant.

:::

Clearly, there are groupings in the data. These groupings could partly be driven by participant demographics (e.g., some participants having larger or smaller regional volumes), but they may also be driven by different structures (e.g., the relationship between surface area and volume may differ across regions). Let's zoom in on two of the larger clusters, and color the points by structure.

:::  {#fig-pairs2}

```{r pairs2}

a2009s |>
  filter(between(SurfArea, 4000, 6000), between(GrayVol, 10000, 20000)) |>
  ggplot(aes(x = SurfArea, y = GrayVol, color = StructName)) +
  geom_point()

```

The data are plotted as in @fig-pairs, but with a restricted axis. Note the apparent variability in the relationship between surface area and volume. 

:::

The relationship between surface area and volume appears to vary by region. For example, the Superior Temporal Sulcus has a similar surface area to the Superior Frontal Gyrus, but that gyrus has a larger volume. 

## Considerations While Working on Projects

### Pivoting the Data

The data have been shared in a ["longer" format](https://r4ds.hadley.nz/data-tidy.html#sec-tidy-data), with rows corresponding to region-level observations. Some analyses benefit from having the data in a ["wider" format](https://r4ds.hadley.nz/data-tidy.html#widening-data), with rows corresponding to participant-level observations. For example, if we're predicting age from cortical thickness, we may want to drop the other morphological measures, and pivot the table such that columns correspond to regional thickness. 

::: {#tbl-wide}

```{r wider}
a2009s |>
  select(StructName, ThickAvg, sub, hemisphere) |>
  pivot_wider(
    names_from = c(StructName, hemisphere),
    values_from = ThickAvg
  ) |>
  head()
```

Average Thickness in a Wider Format. Rows correspond to participants. Note that the region (`StructName`) has been combined with the hemisphere label.

:::

### Variability Across Scanners

{{< include _snippets/mri-scanner-variability.qmd >}}

### Data Quality

{{< include _snippets/mri-qc.qmd >}}

Currently, no QC has been performed on the underlying segmentations and parcellations. With adequate sample sizes, results based on whole-atlas FreeSurfer measurements are generally robust, although there are individual regions that may warrant close inspection (e.g., [McCarthy et al. 2015](http://dx.doi.org/10.3389/fnins.2015.00379), [Vahermaa et al. 2023](https://doi.org/10.1016/j.neuroimage.2023.120306)). Moreover, it remains unclear the extent to which FreeSurfer's accuracy persists over the wide range of demographics found in the A2CPS dataset. The Imaging DIRC is actively exploring these aspects of quality. 

### Methods

FreeSurfer outputs were generated by the [fMRIPrep pipeline](https://fmriprep.org/en/stable/workflows.html#surface-preprocessing). The fMRIPrep pipeline is similar to a typical call to FreeSurfer's `recon-all`, with a majority of changes aiming to parallelize more parts of `recon-all`. However, there are differences, such as in how the brain is masked. The A2CPS use of fMRIPrep is even more different because it entails replacing the fMRIPrep masking procedure with [`SynthStrip`](https://surfer.nmr.mgh.harvard.edu/docs/synthstrip/). In general, these changes are expected to only improve the quality of the FreeSurfer outputs. However, they may hinder some comparisons between results derived from A2CPS and those in other studies.

### Citations

The FreeSurfer package has been developed through extensive research. If you use these derivatives in your analyses, please follow the documentation for citing FreeSurfer: https://surfer.nmr.mgh.harvard.edu/fswiki/FreeSurferMethodsCitation.

{{< include _snippets/a2cps_citations.qmd >}}

When using neuroimaging derivatives, please also cite @sadil_acute_2024.