---
title: Data Structures and Sampling Frames
author: Bradley R. Buchsbaum
date: '`r Sys.Date()`'
output:
  rmarkdown::html_vignette:
    toc: yes
    toc_depth: 2.0
    css: albers.css
    includes:
      in_header: albers-header.html
params:
  family: red
  preset: homage
resource_files:
- albers.css
- albers.js
- albers-header.html

vignette: |
  %\VignetteIndexEntry{Data Structures and Sampling Frames}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r setup, include=FALSE}
if (requireNamespace("ggplot2", quietly = TRUE) && requireNamespace("albersdown", quietly = TRUE)) ggplot2::theme_set(albersdown::theme_albers(family = params$family, preset = params$preset))
knitr::opts_chunk$set(collapse = TRUE, comment = "#>", fig.width=7, fig.height=5, message=FALSE, warning=TRUE)
library(fmrireg)
library(fmridataset) # The fmri_frame container and its accessors
library(neuroim2) # For creating example NeuroVec/NeuroVol
```

```{r albers-classes, echo=FALSE, results='asis'}
cat(sprintf(
  paste0(
    '<script>document.addEventListener("DOMContentLoaded",function(){',
    'document.body.classList.remove("palette-red","palette-lapis","palette-ochre","palette-teal","palette-green","palette-violet","preset-homage","preset-study","preset-structural","preset-adobe","preset-midnight");',
    'document.body.classList.add("palette-%s","preset-%s");',
    '});</script>'
  ),
  params$family,
  params$preset
))
```

## Introduction: Linking Data and Design

Effective fMRI analysis requires associating the measured brain activity (the imaging data) with crucial metadata, including:

*   **Temporal Structure:** When each scan was acquired (TR) and how scans are grouped into runs.
*   **Spatial Structure:** Which brain locations (voxels) are included in the analysis (mask).
*   **Experimental Design:** Timing and properties of experimental events or conditions.

`fmrireg` represents all of this with one container, the `fmri_frame` from the
companion `fmridataset` package. A frame is an *observations by features*
matrix (scans by voxels, vertices, parcels, or components) whose acquisition
timing lives in the observation metadata, whose spatial identity lives in a
typed *feature space*, and whose experimental events live in a keyed event
table. Every modeling function in `fmrireg` (`fmri_lm`, `estimate_betas`,
`estimate_hrf`, `fmri_latent_lm`, ...) takes an `fmri_frame` as its dataset.

This vignette describes how to build a frame from the data you already have and
how to read its timing, events, and spatial layout back. For a first complete
fit, read `vignette("fmrireg", package = "fmrireg")`; return here when you need
to change how scans are stored or divided into runs. For the frame model
itself (views, laziness, persistence) see `vignette("fmridataset", package =
"fmridataset")`.

```{r dataset-decision-table, echo=FALSE}
knitr::kable(data.frame(
  data_you_have = c(
    "NeuroVec runs already in memory",
    "NIfTI runs and a mask on disk",
    "Time-by-feature matrix",
    "Latent components plus spatial loadings",
    "fMRIPrep derivatives (BIDS)"
  ),
  constructor = c(
    "neurovec_frame()",
    "nifti_frame()",
    "matrix_frame()",
    "latent_frame()",
    "fmridataset::read_bids_bold()"
  ),
  feature_space = c(
    "volume_space", "volume_space", "index_space", "basis_space", "volume_space"
  ),
  choose_it_when = c(
    "The full volumes comfortably fit in memory",
    "Runs are large and should be read lazily",
    "Data are ROI, surface, or already tabular",
    "A fmristore LatentNeuroVec already exists",
    "One subject's preprocessed BOLD runs are on disk"
  )
), caption = "Choose a frame constructor from the storage you already have")
```

## The `sampling_frame`

Recall the `sampling_frame` object (introduced in the Overview and detailed in
other vignettes). It defines the temporal structure that the design machinery
expects:

*   `blocklens`: A vector specifying the number of scans (time points) in each run.
*   `TR`: The repetition time (time between scans) in seconds.

```{r sampling_frame_recap}
sframe_example <- sampling_frame(blocklens = c(150, 160), TR = 2.0)
print(sframe_example)
```

A frame does not store a `sampling_frame`. It stores the same facts as two
observation columns, `run_id` and `TR`, and `fmridataset::as_sampling_frame()`
rebuilds the `sampling_frame` from them whenever the design code needs one.
Because the columns are the truth, the timing follows subsetting and
reordering of the frame for free.

## In-Memory Volumetric Data (`neurovec_frame`)

Use this when your fMRI runs are loaded as `neuroim2::NeuroVec` objects in your R session.

**Key Arguments:**

*   `scans`: A `NeuroVec`, or a *list* of them, one for each run.
*   `mask`: A `neuroim2::LogicalNeuroVol` (or a `NeuroVol` whose non-zero voxels form the mask).
*   `TR`: Repetition time (seconds), one value or one per run.
*   `run_length` (Optional): Vector of run lengths; if omitted, inferred from the fourth dimension of each scan.
*   `event_table` (Optional): A `data.frame` containing experimental design information.
*   `censor` (Optional): One entry per scan marking volumes to exclude (logical, 0/1, or 1-based indices).

```{r mem_dataset_example}
# Create minimal example data (2 runs)
d <- c(5, 5, 5, 20) # Small dimensions for example
mask_vol <- neuroim2::LogicalNeuroVol(array(TRUE, d[1:3]), neuroim2::NeuroSpace(d[1:3]))

scan1 <- neuroim2::NeuroVec(array(rnorm(prod(d)), d), neuroim2::NeuroSpace(d))
scan2 <- neuroim2::NeuroVec(array(rnorm(prod(d)), d), neuroim2::NeuroSpace(d))

# Example event table
events_df <- data.frame(
  onset = c(5, 15, 5, 15),
  condition = factor(c("A", "B", "A", "B")),
  run = c(1, 1, 2, 2)
)

# Create the frame
mem_frame <- neurovec_frame(scans = list(scan1, scan2),
                            mask = mask_vol,
                            TR = 2.0,
                            # run_length automatically inferred as c(20, 20)
                            event_table = events_df)

dim(mem_frame)              # scans by masked voxels
class(fmridataset::space(mem_frame))     # a volume_space on the mask's grid
as_sampling_frame(mem_frame)
```

The masked voxel series are extracted once, with `neuroim2::series()`, so a
`SparseNeuroVec` does not need to be densified.

## File-Based Volumetric Data (`nifti_frame`)

This is often the most practical option for typical fMRI analyses where data resides in files.

**Key Arguments:**

*   `scans`: A character vector of file paths to the 4D fMRI image files (e.g., `.nii.gz`), one path per run.
*   `mask`: The path to a 3D mask image (or a `LogicalNeuroVol`, or a `fmridataset::volume_space`).
*   `TR`: Repetition time (seconds).
*   `run_length` (Optional): Volumes per run; if omitted, read from each file's header.
*   `event_table` (Optional): A `data.frame` with experimental design info.
*   `base_path` (Optional): A path to prepend to relative file paths in `scans` and `mask`.

```{r file-dataset-fixture, include=FALSE}
# --- Create Dummy Files (for illustration only) ---
# In a real analysis, these files would already exist.
tmp_dir <- tempdir()
mask_filename <- "mask.nii.gz"
scan1_filename <- "run1.nii.gz"
scan2_filename <- "run2.nii.gz"

mask_file_full_path <- file.path(tmp_dir, mask_filename)
scan1_file_full_path <- file.path(tmp_dir, scan1_filename)
scan2_file_full_path <- file.path(tmp_dir, scan2_filename)

# Create small dummy mask and scans using neuroim2 functionality
d <- c(5, 5, 5) # Mask dimensions
d_run1 <- c(d, 20) # Run 1 dimensions (time=20)
d_run2 <- c(d, 25) # Run 2 dimensions (time=25)

mask_vol_dummy <- neuroim2::NeuroVol(array(1, d), neuroim2::NeuroSpace(d))
scan1_dummy <- neuroim2::NeuroVec(array(rnorm(prod(d_run1)), d_run1), neuroim2::NeuroSpace(d_run1))
scan2_dummy <- neuroim2::NeuroVec(array(rnorm(prod(d_run2)), d_run2), neuroim2::NeuroSpace(d_run2))

# Ensure dummy files are written using their full paths
neuroim2::write_vol(mask_vol_dummy, mask_file_full_path)
neuroim2::write_vec(scan1_dummy, scan1_file_full_path)
neuroim2::write_vec(scan2_dummy, scan2_file_full_path)
# --- End Dummy File Creation ---
```

The hidden fixture above creates tiny valid NIfTI files so the constructor
below genuinely executes. In a study, start with your existing run and mask
paths; the important teaching code is the constructor itself.

```{r file-dataset-example}
# Pass only filenames to 'scans' and 'mask', and specify the directory in 'base_path'
file_frame <- nifti_frame(scans = c(scan1_filename, scan2_filename),
                          mask = mask_filename,
                          TR = 1.5,
                          run_length = c(20, 25), # Must match time dim of files
                          event_table = events_df,
                          base_path = tmp_dir)    # Set base_path to the temp directory

dim(file_frame)
temporal_schema(file_frame)$run_lengths
```

The frame wraps a `fmridataset::nifti_array_source()`: headers and the mask are
read at construction, but no BOLD volumes are read until a fit needs them, and
reads are pushed down per file and per masked voxel. If the files change on
disk after construction, the next read fails with a structured stale-source
error rather than silently returning different values.

```{r file-dataset-cleanup, include=FALSE}
file.remove(mask_file_full_path, scan1_file_full_path, scan2_file_full_path)
```

## Matrix Data (`matrix_frame`)

Use this if your fMRI data is already represented as a 2D matrix where rows are time points and columns are voxels or components (e.g., after surface projection or ROI averaging).

**Key Arguments:**

*   `datamat`: The numeric matrix (time x features).
*   `TR`: Repetition time (seconds).
*   `run_length`: Vector specifying the number of rows (time points) belonging to each run.
*   `event_table` (Optional): One row per experimental event, with onset,
    condition, run, and any modulators used by the model. Its row count usually
    differs from the number of scans; onsets and run labels must instead be
    valid for the declared run lengths.
*   `feature_ids` (Optional): Stable feature IDs; unique column names are used when present.

```{r matrix_dataset_example}
# Example matrix (100 time points, 50 features/voxels)
# Two runs of 50 time points each
time_points <- 100
features <- 50
run_len <- c(50, 50)
example_matrix <- matrix(rnorm(time_points * features), time_points, features)

# Example event table for matrix data
events_mat_df <- data.frame(
  onset = c(seq(5, 45, by=10), seq(5, 45, by=10)),
  condition = factor(rep(c("C", "D"), 10)),
  run = rep(1:2, each = 5)
)

mat_frame <- matrix_frame(datamat = example_matrix,
                          TR = 2.5,
                          run_length = run_len,
                          event_table = events_mat_df)

dim(mat_frame)
class(fmridataset::space(mat_frame))   # an index_space: no spatial geometry
as_sampling_frame(mat_frame)
```

For `matrix_frame`, the concept of a spatial mask is implicit; all columns
provided in `datamat` are included, and the feature space is a plain
`index_space`.

## Latent Data (`latent_frame`)

This constructor is for data that has undergone dimensionality reduction (e.g., PCA, ICA). It takes a `fmristore::LatentNeuroVec`, which stores the basis (latent components over time) and loadings (spatial maps of components), and builds a frame whose assay holds the component scores and whose feature space is a `fmridataset::basis_space` carrying the loadings as its synthesis operator.

**Key Arguments:**

*   `x`: A `LatentNeuroVec` object from the `fmristore` package.
*   `TR`: Repetition time (seconds).
*   `run_length`: Vector specifying run lengths (must sum to the time dimension of `x`).
*   `event_table` (Optional): Experimental design `data.frame`.

```{r latent_dataset_example, eval=FALSE}
# Conceptual example (requires fmristore package and a LatentNeuroVec)
# Assuming 'my_latent_neuro_vec' is a LatentNeuroVec object representing
# 20 components over 300 time points (2 runs of 150)

latent_frame_obj <- latent_frame(
  my_latent_neuro_vec,
  TR = 2.0,
  run_length = c(150, 150),
  event_table = some_event_df
)
collect_assay(latent_frame_obj)              # scores: time x components
basis_synthesis(fmridataset::space(latent_frame_obj))     # loadings: voxels x components
```

Model fitting happens in component space; `fmri_latent_lm()` and the
`latent_sketch` engine reconstruct voxel-wise coefficients through the loadings
when you ask for them.

## Reading a frame back

Once created, a frame is the primary data input for `fmrireg`'s modeling
functions, and the same public accessors work whatever the storage:

```{r accessors}
# Timing: derived from the run_id / TR observation columns
schema <- temporal_schema(mem_frame)
schema$run_lengths
schema$TR
blocklens(as_sampling_frame(mem_frame))

# Events: the keyed event table (an event_id column is added when absent)
head(event_data(mem_frame$tables$events))

# Data: the dense scans-by-features matrix, read under a memory budget
dim(collect_assay(mem_frame))

# Space: the typed feature space and, for volumes, a map back to the grid
n_features(fmridataset::space(mem_frame))
first_volume <- spatial_map(mem_frame, 1)
class(first_volume)
```

Modeling functions take the frame directly:

*   `event_model(..., sampling_frame = as_sampling_frame(frame))`
*   `baseline_model(..., sframe = as_sampling_frame(frame))`
*   `fmri_lm(formula, block, dataset = frame)`
*   `estimate_betas(frame, fixed = ..., ran = ..., block = ...)`

Frames are lazy where they can be. `filter_obs(frame, run_id == "run-2")`
and `frame[, 1:100]` return views that share the source and read nothing
until `collect_assay()` or a fit asks for values.

## Next

- `vignette("a_08_simulation", package = "fmrireg")` — Simulating fMRI data
- `vignette("a_09_linear_model", package = "fmrireg")` — fMRI Linear Model (GLM)
- `vignette("fmridataset", package = "fmridataset")` — Frames, views, and the temporal contract