Four dimensional neuroimaging data

2025-02-17

Working with neuroimaging time-series data

The neuroim2 package contains data structures and functions for reading, accessing, and processing 4-dimensional neuroimaging data.

Reading a four-dimensional NifTI image with read_vec

Here we read in a 4D image consisting of 5 time points,

      library(purrr)
      library(ggplot2)
      file_name <- system.file("extdata", "global_mask_v4.nii", package="neuroim2")
      vec <- read_vec(file_name)
      dim(vec)
#> [1] 64 64 25  4
      vec
#> 
#> DenseNeuroVec (3.13 bytes MB)
#> 
#> - Spatial Info ---------------------------
#> | Dimensions    : 64 x 64 x 25 (4 timepoints)
#> | Total Voxels  : 102,400
#> | Spacing       : 3.5 x 3.5 x 3.7
#> 
#> - Properties ---------------------------
#> | Origin        : 112 x -108 x -46.2
#> | Orientation   : Right-to-Left Posterior-to-Anterior Inferior-to-Superior
#> 
#> - Statistics ---------------------------
#>     Mean +/- SD    : 0.288 +/- 0.453
#> 
#> Label: /home/runner/work/_temp/Library/neuroim2/extdata/global_mask_v4.nii

Now imagine we have a set of 4d images. We can read them in with read_vec. (Here we are just using three versions of the same file for the example).

    
      file_name <- system.file("extdata", "global_mask_v4.nii", package="neuroim2")
      vec <- read_vec(c(file_name, file_name, file_name))
      dim(vec)
#> [1] 64 64 25 12
      
      vec2 <- read_vec(rep(file_name, 10))
      vec2
#> 
#> NeuroVecSeq (10 vectors)
#> 
#> += Sequence Info ---------------------------
#> | Length        : 10
#> | Total Time    : 40 points
#> 
#> += Spatial Info ---------------------------
#> | Dimensions    : 64 x 64 x 25
#> | Spacing       : 3.5 x 3.5 x 3.7
#> | Origin        : 112 x -108 x -46.2
#> | Orientation   : Right-to-Left Posterior-to-Anterior Inferior-to-Superior
#> 
#> += Vector Details --------------------------
#>   1. DenseNeuroVec (4 timepoints)
#>   2. DenseNeuroVec (4 timepoints)
#>   3. DenseNeuroVec (4 timepoints)
#>   4. DenseNeuroVec (4 timepoints)
#>   5. DenseNeuroVec (4 timepoints)
#>   6. DenseNeuroVec (4 timepoints)
#>   7. DenseNeuroVec (4 timepoints)
#>   8. DenseNeuroVec (4 timepoints)
#>   9. DenseNeuroVec (4 timepoints)
#>   10. DenseNeuroVec (4 timepoints)

To extract a subset of volumes we can use the sub_vector function:

    
      vec_1_6 <- sub_vector(vec, 1:6)
      dim(vec_1_6)
#> [1] 64 64 25  6
      vec_1_6
#> 
#> NeuroVecSeq (2 vectors)
#> 
#> += Sequence Info ---------------------------
#> | Length        : 2
#> | Total Time    : 6 points
#> 
#> += Spatial Info ---------------------------
#> | Dimensions    : 64 x 64 x 25
#> | Spacing       : 3.5 x 3.5 x 3.7
#> | Origin        : 112 x -108 x -46.2
#> | Orientation   : Right-to-Left Posterior-to-Anterior Inferior-to-Superior
#> 
#> += Vector Details --------------------------
#>   1. DenseNeuroVec (4 timepoints)
#>   2. DenseNeuroVec (2 timepoints)

Extracting time-series data using the series and series_roi functions

To get the time-series at voxel (1,1,1) we can use the series function:

      
      series(vec_1_6, 1,1,1)
#> [1] 0 0 0 0 0 0

We can extract a 4d region of interest with the series_roi as follows:

      file_name <- system.file("extdata", "global_mask_v4.nii", package="neuroim2")
      vol <- read_vol(file_name)
      roi <- spherical_roi(vol, c(12,12,12), radius=8)
      rvec1 <- series_roi(vec, roi)
      
      ## or alternatively as a pipeline
      rvec2 <- read_vol(file_name) %>% spherical_roi(c(12,12,12), radius=8) %>% series_roi(vec,.)
      rvec2
#> 
#>  === ROIVec Object === 
#> 
#> - Structure 
#>   Points:     49
#>   Features:   3 (147 total)
#>   Memory:     12.7 Kb
#> 
#> - Spatial Properties
#>   Parent Space: 64 x 64 x 25 x 12
#>   Centroid:     [13.0, 13.0, 13.0 mm]
#> 
#> - Value Properties
#>   Range:    [0.00, 0.00]
#> 
#> ======================================
#> 
#>  Access Methods: 
#>   .  Get Points:   coords(object) 
#>   .  Get Values:   as.matrix(object) 
#>   .  Subset:       object[1:10, ]
      
      ## we can extract the ROI values with the `values` method.
      assertthat::assert_that(all(values(rvec1) == values(rvec2)))
#> [1] TRUE
      assertthat::assert_that(all(coords(rvec1) == coords(rvec2)))
#> [1] TRUE

We can also extract an ROI using 1d indices:

r1 <- series_roi(vec, 1:100)
r1
#> 
#>  === ROIVec Object === 
#> 
#> - Structure 
#>   Points:     100
#>   Features:   3 (300 total)
#>   Memory:     17.5 Kb
#> 
#> - Spatial Properties
#>   Parent Space: 64 x 64 x 25 x 12
#>   Centroid:     [27.5, 1.4, 1.0 mm]
#> 
#> - Value Properties
#>   Range:    [0.00, 0.00]
#> 
#> ======================================
#> 
#>  Access Methods: 
#>   .  Get Points:   coords(object) 
#>   .  Get Values:   as.matrix(object) 
#>   .  Subset:       object[1:10, ]

Or we can extract a plain matrix using the series function:

r2 <- series(vec, 1:100)
dim(r2)
#> [1]  12 100

We can also use coordinate indexing using voxel coordinates. First we load a binary mask with the same spatial dimensions as our NeuroVec:

mask <- read_vol(system.file("extdata", "global_mask_v4.nii", package="neuroim2"))

Now we convert indices to voxels and extract a matrix of values at the specified locations:

vox <- index_to_grid(mask, 1:100)

r3 <- series(vec, vox)
dim(r3)
#> [1]  12 100

And the same using series_roi:

r4 <- series_roi(vec,vox)
r4
#> 
#>  === ROIVec Object === 
#> 
#> - Structure 
#>   Points:     100
#>   Features:   3 (300 total)
#>   Memory:     18.7 Kb
#> 
#> - Spatial Properties
#>   Parent Space: 64 x 64 x 25 x 12
#>   Centroid:     [27.5, 1.4, 1.0 mm]
#> 
#> - Value Properties
#>   Range:    [0.00, 0.00]
#> 
#> ======================================
#> 
#>  Access Methods: 
#>   .  Get Points:   coords(object) 
#>   .  Get Values:   as.matrix(object) 
#>   .  Subset:       object[1:10, ]