Spatially Regularized Super-Resolved Constrained Spherical Deconvolution (SR$^2$-CSD) of Diffusion MRI Data

TL;DR

SR$^2$-CSD introduces a spatially regularized, super-resolved approach to constrained spherical deconvolution by embedding a TV-denoised, J-invariance-calibrated prior into a nonnegativity-constrained quadratic program. This framework enables higher angular resolution (e.g., $l_{max}=12$) while improving stability, spatial coherence, and reproducibility of FOD estimates across phantoms, in vivo data, and simulated connectivity, with higher fidelity to ground-truth connectivity in tractography. Across diverse datasets, SR$^2$-CSD consistently lowers angular and peak errors, enhances repeatability, and yields connectivity metrics more aligned with ground truth, often with statistical significance after correction. The work demonstrates the feasibility and practical benefits of incorporating spatial priors into CSD, and provides an open-source implementation in DIPY for broader diffusion MRI research and clinical applicability.

Abstract

Constrained Spherical Deconvolution (CSD) is widely used to estimate the white matter fiber orientation distribution (FOD) from diffusion MRI data. Its angular resolution depends on the maximum spherical harmonic order ($l_{max}$): low $l_{max}$ yields smooth but poorly resolved FODs, while high $l_{max}$, as in Super-CSD, enables resolving fiber crossings with small inter-fiber angles but increases sensitivity to noise. In this proof-of-concept study, we introduce Spatially Regularized Super-Resolved CSD (SR$^2$-CSD), a novel method that regularizes Super-CSD using a spatial FOD prior estimated via a self-calibrated total variation denoiser. We evaluated SR$^2$-CSD against CSD and Super-CSD across four datasets: (i) the HARDI-2013 challenge numerical phantom, assessing angular and peak number errors across multiple signal-to-noise ratio (SNR) levels and CSD variants (single-/multi-shell, single-/multi-tissue); (ii) the Sherbrooke in vivo dataset, evaluating spatial coherence of FODs; (iii) a six-subject test-retest dataset acquired with both full (96 gradient directions) and subsampled (45 directions) protocols, assessing reproducibility; and (iv) the DiSCo phantom, evaluating tractography accuracy under varying SNR levels and multiple noise repetitions. Across all evaluations, SR$^2$-CSD consistently reduced angular and peak number errors, improved spatial coherence, enhanced test-retest reproducibility, and yielded connectivity matrices more strongly correlated with ground-truth. Most improvements were statistically significant under multiple-comparison correction. These results demonstrate that incorporating spatial priors into CSD is feasible, mitigates estimation instability, and improves FOD reconstruction accuracy.

Figures (15)

• Figure 1: Selection of the TV denoising strength $K$ for SR$^2$-CSD with $l_{max} = 12$ in the HARDI Phantom (SNR = 15). (a) Angular Error (AE). (b) Peak Number Error (PNE). Horizontal lines show baseline performance of standard CSD ($l_{max} = 8$, blue) and Super-CSD ($l_{max} = 12$, orange). The vertical red line indicates the $K$ value automatically selected via J-invariance, which matches the minimum of the curves.
• Figure 2: Reconstruction accuracy on the HARDI phantom (SNR = 30, 10 noise realizations, for single-shell single-tissue reconstruction). (a) Without MPPCA denoising. (b) After MPPCA denoising. Boxplots show Angular Error (AE) and Peak Number Error (PNE) for CSD ($l_{max}=8$), Super-CSD ($l_{max}=12$), and SR$^2$-CSD ($l_{max}=12$). Statistical significance assessed via paired t-tests with Bonferroni correction. Asterisks indicate significance levels: *$p < 0.05$; **$p < 0.01$; ***$p < 0.001$.
• Figure 3: Reconstruction accuracy on the HARDI phantom (SNR = 30, 10 noise realizations, for multi-shell multi-tissue reconstruction). (a) Without MPPCA denoising. (b) After MPPCA denoising. Boxplots show Angular Error (AE) and Peak Number Error (PNE) for CSD ($l_{max}=8$), Super-CSD ($l_{max}=12$), and SR$^2$-CSD ($l_{max}=12$). Statistical significance assessed via paired t-tests with Bonferroni correction. Asterisks indicate significance levels: *$p < 0.05$; **$p < 0.01$; ***$p < 0.001$.
• Figure 4: Visual comparison of estimated Fiber Orientation Distributions (FODs) in two regions of the HARDI Phantom with fiber crossings. Data: single-shell ($b=3000;s/\text{mm}^2$), SNR = 15, after MPPCA denoising. Row 1: ground-truth fiber orientations. Rows 2–3: estimated FODs from CSD ($l_{max}=8$), Super-CSD ($l_{max}=12$), and SR$^2$-CSD ($l_{max}=12$). Three regions are highlighted where SR$^2$-CSD improved the angular resolution (red boxes) and the angular coherency (blue box).
• Figure 5: Fiber Orientation Distributions (FODs) estimated from the raw Sherbrooke data and MPPCA-denoised data acquired with a single shell ($b = 3500$ s/mm$^2$, 64 gradient directions). From left to right, the figure shows results for CSD with $l_{\text{max}} = 8$, Super-CSD with $l_{\text{max}} = 12$, and SR$^2$-CSD with $l_{\text{max}} = 12$. The red boxes highlight regions where denoising improves the FOD estimates from CSD and Super-CSD, yielding more spatially coherent lobes and reducing isolated peaks that do not appear in neighboring voxels. The blue box marks an adjacent superior region where, even after denoising, the anterior–posterior green lobes visible in the red boxes become markedly attenuated or nearly disappear. In contrast, SR$^2$-CSD preserves these lobes across both regions, especially in the denoised data, producing greater spatial continuity.