Smoothing and filtering with neuroim2

This vignette gives a practical overview of the main spatial and spatio‑temporal smoothing tools in neuroim2:

We’ll work with the demo volume used elsewhere in the package. All code assumes this setup:

set.seed(1)

demo_path <- system.file("extdata", "global_mask_v4.nii", package = "neuroim2")
vol3d <- read_vol(demo_path)          # 3D volume
vec4d <- read_vec(demo_path)          # 4D NeuroVec (same file)

sp3  <- space(vol3d)
dims <- dim(vol3d)

1. Gaussian smoothing with gaussian_blur()

gaussian_blur() applies an isotropic Gaussian kernel within an optional mask. It is useful as a simple, fast baseline smoother.

blur_light <- gaussian_blur(vol3d, vol3d, sigma = 2, window = 1)
blur_strong <- gaussian_blur(vol3d, vol3d, sigma = 4, window = 2)

dim(blur_light)
#> [1] 64 64 25

Use smaller sigma and window for gentle smoothing; larger values blur more but can wash out small structures.

2. Edge‑preserving guided filter with guided_filter()

guided_filter() smooths while preserving edges by fitting local linear models between the input and output.

gf_vol <- guided_filter(vol3d, radius = 4, epsilon = 0.7^2)
gf_vol
#> 
#>  === NeuroVol Object === 
#> 
#> * Basic Information 
#>   Type:      DenseNeuroVol
#>   Dimensions: 64 x 64 x 25 (806.6 Kb)
#>   Total Voxels: 102,400
#> 
#> * Data Properties
#>   Value Range: [0.00, 1.00]
#> 
#> * Spatial Properties
#>   Spacing: 3.50 x 3.50 x 3.70 mm
#>   Origin:  112.0, -108.0, -46.2 mm
#>   Axes:    Right-to-Left x Posterior-to-Anterior x Inferior-to-Superior
#> 
#> ======================================
#> 
#>  Access Methods: 
#>   .  Get Slice:   slice(object, zlevel=10) 
#>   .  Get Value:   object[i, j, k] 
#>   .  Plot:       plot(object)  # shows multiple slices

Compared to Gaussian smoothing, guided filtering tends to retain sharp boundaries between tissues while denoising within regions.

3. Bilateral spatial filter with bilateral_filter()

bilateral_filter() combines spatial distance and intensity similarity, reducing noise while respecting edges.

bf_vol <- bilateral_filter(
  vol3d,
  spatial_sigma   = 2,
  intensity_sigma = 1,
  window          = 1
)
bf_vol
#> 
#>  === NeuroVol Object === 
#> 
#> * Basic Information 
#>   Type:      DenseNeuroVol
#>   Dimensions: 64 x 64 x 25 (806.6 Kb)
#>   Total Voxels: 102,400
#> 
#> * Data Properties
#>   Value Range: [0.00, 0.99]
#> 
#> * Spatial Properties
#>   Spacing: 3.50 x 3.50 x 3.70 mm
#>   Origin:  112.0, -108.0, -46.2 mm
#>   Axes:    Right-to-Left x Posterior-to-Anterior x Inferior-to-Superior
#> 
#> ======================================
#> 
#>  Access Methods: 
#>   .  Get Slice:   slice(object, zlevel=10) 
#>   .  Get Value:   object[i, j, k] 
#>   .  Plot:       plot(object)  # shows multiple slices

Key parameters:

Smaller intensity_sigma better preserves edges; larger values behave more like pure Gaussian smoothing.

4. Laplacian enhancement with laplace_enhance()

laplace_enhance() is designed for sharpening rather than smoothing. It uses a multi‑layer, non‑local Laplacian scheme to enhance details.

sharp_vol <- laplace_enhance(vol3d, k = 2, patch_size = 3,
                             search_radius = 2, h = 0.7)
sharp_vol
#> 
#>  === NeuroVol Object === 
#> 
#> * Basic Information 
#>   Type:      DenseNeuroVol
#>   Dimensions: 64 x 64 x 25 (806.6 Kb)
#>   Total Voxels: 102,400
#> 
#> * Data Properties
#>   Value Range: [-0.02, 1.01]
#> 
#> * Spatial Properties
#>   Spacing: 3.50 x 3.50 x 3.70 mm
#>   Origin:  112.0, -108.0, -46.2 mm
#>   Axes:    Right-to-Left x Posterior-to-Anterior x Inferior-to-Superior
#> 
#> ======================================
#> 
#>  Access Methods: 
#>   .  Get Slice:   slice(object, zlevel=10) 
#>   .  Get Value:   object[i, j, k] 
#>   .  Plot:       plot(object)  # shows multiple slices

Use higher h or more layers k for stronger enhancement, but be cautious about amplifying noise.

5. 4D bilateral smoothing with bilateral_filter_4d()

For 4D time‑series data (NeuroVec), bilateral_filter_4d() extends the bilateral idea across space and time:

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

bf4d <- bilateral_filter_4d(
  vec4d, mask3d,
  spatial_sigma   = 2,
  intensity_sigma = 1,
  temporal_sigma  = 1,
  spatial_window  = 1,
  temporal_window = 1,
  temporal_spacing = 1
)

dim(bf4d)
#> [1] 64 64 25  4

This is useful when you want joint spatial + temporal denoising that still respects intensity boundaries (e.g. fMRI or 4D structural sequences).

6. Correlation‑guided bilateral filtering with cgb_filter()

cgb_filter() implements correlation‑guided bilateral (CGB) smoothing: it builds a sparse graph based on spatial distance and time‑series correlations, then diffuses data over that graph.

Basic usage mirrors the bilateral interface:

cgf <- cgb_filter(
  vec4d,
  mask          = mask3d,
  spatial_sigma = 3,
  window        = NULL,      # auto from spatial_sigma and spacing
  corr_map      = "power",
  corr_param    = 2,
  topk          = 16,
  passes        = 1,
  lambda        = 1
)

dim(cgf)
#> [1] 64 64 25  4

Parameter intuition:

Internally this calls cgb_make_graph() to build the graph and then cgb_smooth() to diffuse over it. Use return_graph = TRUE to keep the graph for reuse.

cg_out <- cgb_filter(vec4d, mask3d,
                     spatial_sigma = 3, window = NULL,
                     topk = 16, return_graph = TRUE)

str(cg_out$graph)
#> List of 5
#>  $ row_ptr : int [1:29533] 0 1 2 3 4 5 6 7 8 9 ...
#>  $ col_ind : int [1:29532] 474 475 476 477 480 481 482 483 536 537 ...
#>  $ val     : num [1:29532] 1 1 1 1 1 1 1 1 1 1 ...
#>  $ dims3d  : int [1:3] 64 64 25
#>  $ mask_idx: int [1:29532] 475 476 477 478 481 482 483 484 537 538 ...

Choosing a smoother