LibppRPA: Data-Adaptive Brain Parcellation

• LibppRPA is an open-source software package that performs principal parcellation analysis by clustering tractography endpoints to derive interpretable, data-driven brain connectomes.
• It aggregates fiber endpoints from diffusion imaging and applies mini-batch K-means with a bidirectional distance metric to form population-level fiber-bundle bases.
• The package integrates seamlessly with tractography pipelines and sparse regression models to enable efficient, reproducible trait prediction in neuroimaging studies.

LibppRPA is an open-source software package for Principal Parcellation Analysis (PPA), designed to perform tractography-based, data-adaptive parcellation of brain structural connectomes and enable trait-predictive modeling across populations. By moving away from the conventional reliance on atlas-based regions of interest (ROIs), LibppRPA employs clustering of fiber endpoints to generate population-level fiber-bundle bases, yielding lower-dimensional, interpretable compositional representations that facilitate statistical analyses and prediction tasks in neuroimaging (Liu et al., 2021).

1. Mathematical and Statistical Foundations

LibppRPA formalizes brain connectome representation through a sequence of operations:

• Fiber Aggregation: For $n$ subjects, each with $m_i$ reconstructed fibers (via tractography), endpoints are denoted $a_{ik}, b_{ik}\in\mathbb{R}^3$. The aggregate data matrix $Z\in\mathbb{R}^{6\times M}$, with $M=\sum_i m_i$, column-stacks all ordered endpoint pairs across the cohort.
• Clustering (K-means with bidirectional distance): $Z$ columns are clustered using mini-batch $K$-means to obtain $K$ clusters $\mathcal{A}_K = \{A_K^{(1)}, ..., A_K^{(K)}\}$, with centers $c_1, ..., c_K$ solving

$$\min_{c_1,\dots,c_K} \sum_{j=1}^M \min_{1\leq k\leq K} d(z_j,c_k)^2$$

where $$d((a,b),c) = \min(\|[a;b] - c\|_2, \|[b;a] - c\|_2 )$$ ensures invariance to endpoint order.

• Compositional Connectome Encoding: Each subject $i$’s connectome becomes a $K$-vector

$$\omega_{i} = (\omega_{i1}, ..., \omega_{iK}), \qquad \omega_{ik} = \frac{|\{f_{ij} : (a_{ij}, b_{ij}) \in \text{cluster } k\}|}{m_i}$$

resulting in the matrix $\Omega\in\mathbb{R}^{n\times K}$ for downstream analysis.

• Trait Prediction via Sparse Modeling: To relate $\Omega$ to a scalar trait $y_i$, regularized linear models (e.g., LASSO) are fit:

$$\min_{\beta_0, \beta}\ \sum_{i=1}^n (y_i - \beta_0 - \sum_{k=1}^{K-1} \omega_{ik}\beta_k)^2 + \lambda\|\beta\|_1$$

Enforcing $\sum_k\omega_{ik}=1$ (compositional normalization) removes overparameterization.

This approach reduces reliance on a priori atlas definitions and leverages data-derived bundles to capture subject-level connectome structure in a low-dimensional but descriptive manner (Liu et al., 2021).

2. Algorithmic Pipeline

LibppRPA operationalizes the above statistical framework in a unified workflow comprising three modules:

• Module (i): Fiber Reconstruction
  • Inputs: raw diffusion-weighted imaging (DWI) and T1 MRI per subject.
  • Processing: TractoFlow pipeline (Nextflow + Singularity) reconstructs tractograms, typically producing 2–3 million streamlines per subject. Optional outlier removal via QuickBundles is supported.
• Module (ii): Data-adaptive Parcellation
  • Extract all ($a$, $b$) pairs for each fiber.
  • Aggregate endpoint pairs for all subjects into the $6\times M$ matrix $Z$.
  • Apply mini-batch KMeans ($K$ clusters, default batch-size $B=1000$), using the bidirectional distance for fiber symmetry. Clusters define combinatorial fiber-bundle “parcels.”
  • For each subject, compute $\omega_i$ via normalized cluster counts of their streamlines.
• Module (iii): Trait-adaptive Supervised Learning
  • Given $(\Omega, y)$, fit regularized regression models such as LASSO or ElasticNet. Nonzero coefficients $\beta_k$ indicate bundles predictive of the trait.

This end-to-end pipeline is optimized for compatibility with high-throughput imaging and machine learning ecosystems, with explicit support for parameter tuning and reproducible batch processing (Liu et al., 2021).

3. Library API and Usage

LibppRPA is implemented as a pure Python package, installable via Conda or pip, but requiring system-level dependencies for tractography (Nextflow, Singularity/Docker, MRtrix3/FSL).

Key API features:

• PPA class encapsulates the workflow:
  • .transform(): returns $\Omega$ ($n\times K$ compositional matrix).
  • .get_cluster_centers(): returns cluster centers $c_k\in\mathbb{R}^6$.
  • .get_assignments(): per-subject array of cluster labels for their streamlines.

Statistical modeling is decoupled and uses standard libraries (scikit-learn/ElasticNet/LASSO) downstream:

from sklearn.linear_model import LassoCV lasso = LassoCV(cv=5).fit(omega, y) active_idx = np.where(lasso.coef_ != 0)[0]

Parameter recommendations:

• $K$: typical range $50$–$500$; select via 5-fold CV on downstream MSE.
• batch_size=1000 balances memory and speed.
• random_state controls reproducibility.

The API structure affords rapid integration into connectome-analysis and trait-mapping pipelines (Liu et al., 2021).

4. Performance Characteristics and Example Analyses

Empirical evaluation on Human Connectome Project (HCP) data ($n=1065$) shows:

• For $K=400$, 5-fold CV mean squared error (MSE) for predicting PicVocab scores was $30$–$40$ points lower than classical APA-based (atlas parcellation analysis) methods (SBL, MultiGraphPCA).
• Model parsimony: PPA typically selects $10$–$20$ nonzero $\beta_k$ for trait prediction, versus hundreds in APA approaches.
• Consistency: Performance is robust to tractography algorithm (TractoFlow, EuDX, SFM) and regularization scheme (LASSO, ElasticNet).
• Cross-validated hyperparameter selection over $K$ displays a characteristic U-shaped MSE curve.

Hyperparameter cross-validation example:

def cv_mse_over_K(K_values, dwi_list, bval_list, bvec_list, y): mse_list = [] for K in K_values: ppa = PPA(n_clusters=K, batch_size=1000, ... ) omega = ppa.fit_transform(dwi_list, bval_list, bvec_list) # 5-fold split and LASSO here... mse_list.append(mean_mse) return mse_list Ks = [10, 50, 100, 200, 300, 400, 500] mse_vals = cv_mse_over_K(Ks, dwi_list, bval_list, bvec_list, y)

This demonstrates objective, reproducible performance quantification and model selection (Liu et al., 2021).

5. Extensions, Integration, and Practical Advice

LibppRPA is portable: its intermediate output $\Omega$ is a standard matrix suitable for input into any Python-based machine learning algorithm. Visualization and further analysis are enabled by exporting cluster centers and streamlines for anatomical mapping (DSI Studio, nibabel).

Extensions are possible by:

• Replacing KMeans with alternative clustering (spectral clustering, NMF).
• Utilizing different regression models (ElasticNet, kernel machines, SCAD).
• Visualizing or exporting “active” bundles as discovered by nonzero $\beta_k$.

Practical notes:

• The dependence on tractography toolchains (TractoFlow, QuickBundles, etc.) requires suitable infrastructure (Linux/Mac, containerization).
• Choice of $K$ is critical; CV-guided selection is recommended.
• For visualization or post hoc anatomical interpretation, export cluster centers to .trk/.tck or similar formats.

This modularity enables deep integration with neuroimaging pipelines, facilitating advanced connectome-based analyses (Liu et al., 2021).

6. Impact and Methodological Significance

By breaking dependence on arbitrary ROI atlases and adjacency matrix-based features (which scale as $O(p^2)$ in number of atlas regions), LibppRPA reduces dimensionality to $O(K)$, thereby improving interpretability, statistical power, and computational tractability.

Compared to prior approaches in connectomics:

• The tractography-driven, clustering-based parcellation yields representations adaptable to the population and the trait of interest, rather than imposing brain partitions a priori.
• The compositional encoding is interpretable as fiber-bundle proportions, and the resulting basis can be visualized anatomically.
• Empirical results, including those from HCP, demonstrate strong predictive power and model parsimony for behavioral traits.

A plausible implication is that this approach offers a statistically robust, trait-adaptive means for connectome dimensionality reduction in large-scale neuroimaging studies, and provides a bridge between population-level imaging and hypothesis-driven neuroanatomical inference (Liu et al., 2021).