Recursive Nearest Agglomeration (ReNA) Clustering

Performs spatially-constrained hierarchical clustering using recursive 1-nearest neighbor aggregation. ReNA is a fast, linear-time agglomerative algorithm that avoids the "percolation" problem of standard hierarchical methods by using 1-NN graphs to ensure balanced cluster sizes.

Usage

rena(
  bvec,
  mask,
  K = 100,
  connectivity = 26,
  max_iterations = 50,
  verbose = FALSE,
  exact_k = TRUE
)

Arguments

bvec: A NeuroVec instance supplying the 4D data to cluster.
mask: A NeuroVol mask defining the voxels to include in clustering. Exactly finite values greater than zero are included.
K: The number of clusters to find (default 100).
connectivity: Exact masked-grid connectivity (6, 18, or 26).
max_iterations: Maximum number of recursion iterations (default 50). Algorithm stops when K clusters are reached or max_iterations is hit.
verbose: Logical; whether to print progress messages. Default FALSE.
exact_k: Logical; whether to use the shared adjacency-preserving merge and split engine to ensure exactly K clusters. A target below the number of disconnected mask components is infeasible.

Value

A list of class rena_cluster_result (inheriting from cluster_result) with the following elements:

clusvol: An instance of type ClusteredNeuroVol.
cluster: A vector of cluster indices equal to the number of voxels in the mask.
centers: A matrix of cluster centers with each row representing the feature vector for a cluster.
coord_centers: A matrix of spatial coordinates with each row corresponding to a cluster.
n_clusters: Actual number of clusters found.
method: Clustering method used ("rena").
parameters: List of all parameters used.
metadata: List containing iteration count and convergence information.

Details

Algorithm Description

ReNA recursively aggregates features using a 1-Nearest Neighbor (1-NN) graph to ensure clusters remain roughly balanced in size and spatially compact. The algorithm proceeds in iterations:

Connectivity Constraint: Clustering respects the spatial topology graph G. Distance is only calculated between spatially adjacent voxels.
Similarity Calculation: For all connected edges, compute squared Euclidean distance in feature space.
1-NN Subgraph: For every voxel, identify its single nearest neighbor. This forms a directed subgraph Q. The 1-NN graph is less likely to "percolate" (create giant components) compared to k-NN graphs where k >= 2.
Connected Components: Extract weakly connected components of Q. These components become the clusters for the current iteration.
Reduction: Features within each component are averaged to form a single "super-feature" for the next iteration. The graph G is contracted so edges between new clusters exist if any constituent voxels were connected.
Stopping Condition: Repeat recursively until the number of clusters drops to the desired K. If exact_k = TRUE, edges are pruned by distance to ensure exactly K clusters.

Performance Characteristics

Complexity: O(N log N) per iteration where N = number of voxels
Memory: O(N) for graph and assignments
Speed: Generally faster than iterative methods (supervoxels)
Iterations: Typically 5-15 iterations to reach K clusters from N voxels

Advantages over Other Methods

No percolation: Unlike single-linkage, doesn't form one giant cluster
Balanced sizes: 1-NN graph encourages roughly equal cluster sizes
Topology-aware: Respects spatial structure via connectivity graph
Deterministic: No random initialization
Fast convergence: Linear-time operations per iteration

Parallelization Status

Currently NOT parallelized. ReNA uses a sequential recursive algorithm that processes the graph in a hierarchical manner. Each iteration depends on the previous iteration's clustering result.

Sequential Operations:

Distance Computation: Pairwise distances for connected voxels (C++)
1-NN Finding: Each voxel finds its nearest neighbor (C++)
Component Detection: Union-Find algorithm for connected components (C++)
Feature Aggregation: Mean pooling within components (C++)
Graph Contraction: Building reduced graph for next iteration (C++)

Why Not Parallelized:

Sequential dependency: Each iteration requires previous iteration's result
Fast per-iteration: Linear-time operations make parallelization overhead high
Small iteration count: Typically 5-15 iterations total
C++ optimization: Core operations already optimized in C++

Performance Tips:

Reduce connectivity: Use connectivity=6 instead of 26 for larger graphs
Use smaller K: Fewer target clusters means fewer iterations
Pre-smooth data: Reduces noise, improves cluster coherence
Parallel alternatives: In the unified API, supervoxels and slic support parallel = TRUE. The other methods do not advertise a parallel clustering contract.

Note

Consider using cluster4d with method = "rena" for a standardized interface across all clustering methods.

References

Hoyos-Idrobo, A., Varoquaux, G., Kahn, J., & Thirion, B. (2019). Recursive Nearest Agglomeration (ReNA): Fast clustering for approximation of structured signals. Pattern Recognition, 94, 17-28.

Examples

if (FALSE) { # \dontrun{
# Small synthetic example
library(neuroim2)
mask <- NeuroVol(array(1, c(10,10,10)), NeuroSpace(c(10,10,10)))
vec <- replicate(20,
                 NeuroVol(array(runif(10*10*10), c(10,10,10)),
                          NeuroSpace(c(10,10,10))),
                 simplify=FALSE)
vec <- do.call(concat, vec)

# Run ReNA clustering
rena_res <- rena(vec, mask, K=50, connectivity=26)
print(rena_res$n_clusters)

# With verbose output
rena_res <- rena(vec, mask, K=50, verbose=TRUE)
} # }