# Benchmarking 3D multi-coil NC-PDNET MRI reconstruction

Asma Tanabene, Chaithya Giliyar Radhakrishna, Aurélien Massire, Mariappan S. Nadar, Philippe Ciuciu

<figure><img src="image_1.png" /></figure>

## To cite this version:

Asma Tanabene, Chaithya Giliyar Radhakrishna, Aurélien Massire, Mariappan S. Nadar, Philippe Ciuciu. Benchmarking 3D multi-coil NC-PDNET MRI reconstruction. ISMRM & ISMRT 2025 - Annual Meeting & Exhibition, May 2025, Honolulu, Hawaii, United States. hal-04771143

HAL Id: hal-04771143

https://inria.hal.science/hal-04771143v1

Submitted on 7 Nov 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# BENCHMARKING 3D MULTI-COIL NC-PDNET MRI RECONSTRUCTION

A. Tanabene 1,2,3, Chaithya G. R. 2,3, A. Massire 1, M. Nadar 4, P. Ciuciu 2,3.

¹Siemens Healthineers, Courbevoie, 92400, France

²CEA, NeuroSpin, CNRS, Université Paris-Saclay, Gif-sur-Yvette, 91191, France

³Inria, MIND, Palaiseau, 91120, France

⁴Siemens Healthineers, Princeton, 08540, NJ, USA

## ABSTRACT

Deep learning has shown great promise for MRI reconstruction from undersampled data, yet there is a lack of research on validating its performance in 3D parallel imaging acquisitions with non-Cartesian undersampling. In addition, the artifacts and the resulting image quality depend on the under-sampling pattern. To address this uncharted territory, we extend the Non-Cartesian Primal-Dual Network (NC-PDNet), a state-of-the-art unrolled neural network, to a 3D multi-coil setting. We evaluated the impact of channel-specific versus channel-agnostic training configurations and examined the effect of coil compression. Finally, we benchmark four distinct non-Cartesian undersampling patterns, with an acceleration factor of six, using the publicly available Calgary-Campinas dataset. Our results show that NC-PDNet trained on compressed data with varying input channel numbers achieves an average PSNR of 42.98dB for 1 mm isotropic 32 channel whole-brain 3D reconstruction. With an inference time of 4.95sec and a GPU memory usage of 5.49 GB, our approach demonstrates significant potential for clinical research application.

**Index Terms**— MRI, image reconstruction, deep learning, non-Cartesian, multi-coil 3D imaging.

## 1. INTRODUCTION

Magnetic Resonance Imaging (MRI) is a non-invasive medical imaging technique that offers excellent soft tissue contrast. Advances in MRI technology are constantly evolving to achieve higher image resolution and acquisition speed, enabling faster, more detailed, and more accurate imaging. In the quest for high SNR and high isotropic spatial resolution images, 3D parallel-imaging protocols are increasingly used, often resulting in longer scan times. Inspired by compressed sensing theories [1], accelerated scans can be performed by subsampling the k-space in a variable density. The fully 3D non-Cartesian imaging paradigm allows us to attain such variable density sampling, with denser sampling at the center of the k-space and tapering off towards the periphery, offering faster and more efficient coverage of the k-space than Cartesian sampling [2]. However, while acquisition speed increases, this change also transfers computation complexity and time requirements to the image reconstruction task.

Recently, deep learning has emerged as an alternative for fast MR image reconstruction without compromising image quality. Among the most promising solutions in this field are physics-driven neural networks that unroll optimization algorithms [3, 4, 5, 6] and produce enhanced MR images from raw k-space data. However, most recent deep neural networks do not effectively scale to

3D multi-coil acquisition setups. For example, the non-Cartesian primal-dual network (NC-PDNet) [7], the first density-compensated unrolled neural network designed for non-Cartesian imaging, has only been validated on multi-coil 2D and single-coil magnitude-only 3D k-space data.

Training end-to-end unrolled neural networks on 3D multi-coil non-Cartesian data presents significant challenges. One major issue is that handling non-Cartesian data requires using the Non-Uniform Fast Fourier Transform (NUFFT) operator, which is more computationally and memory-intensive than the standard Fast Fourier Transform, especially as the number of channels increases. Furthermore, incorporating 3D convolutional neural networks (CNN) into unrolled architectures requires managing a much larger number of parameters (i.e. neurons) than in 2D configurations. These high memory requirements have significantly constrained research efforts in scaling recent advances to 3D multi-coil non-Cartesian settings. To manage this, memory-efficient solutions, such as coil compression, are recommended to reduce the input size of 3D multi-coil k-space data.

In this work, we address the challenge of extending NC-PDNet to non-Cartesian 3D multi-coil MRI reconstruction. Further, to understand the artifacts and reconstructed image quality, we conduct a retrospective ablation study with four different non-Cartesian undersampling patterns. Finally, we evaluate the impact of coil compression and compare channel-specific versus channel-agnostic training approaches.

## 2. EXPERIMENTAL SETUP

### 2.1. Model

In our multi-coil non-Cartesian under-sampling context, the reconstructed image is recovered by solving the following ill-posed optimization problem:

$$ \arg\min_{\mathbf{x} \in \mathbb{C}^N} \frac{1}{2L} \sum_{\ell=1}^{L} \|\mathbf{y}_{\ell} - \mathcal{F}_{\Omega} \mathbf{S}_{\ell} \mathbf{x}\|_{2}^{2} + \mathcal{R}(\mathbf{x}) \quad (1) $$

where $\mathcal{F}$ is the multi-coil 3D NUFFT with locations defined by sampling pattern $\Omega$ applied on $L$ number of coils. $\mathbf{S}_{\ell}$ refers to the sensitivity map of the $\ell^{th}$ corresponding coil, and $\mathcal{R}(\mathbf{x})$ serves as a prior for the solution, which is conventionally sparsity enforced in wavelet domain. $\mathbf{y}_{\ell} \in \mathbb{C}^M$ are the k-space measurements and $\mathbf{x} \in \mathbb{C}^N$ is the target image, where $M \ll N$. NC-PDNet unrolls the proximal gradient descent, an iterative algorithm that solves our optimization problem (1), as follows:
$$ 
\begin{aligned}
x_{n+1} &= x_n - \epsilon_n A^H (Ax_n - y), \\
x_{n+1} &= \text{prox}_{\epsilon_n \mathcal{R}} (x_{n+1})
\end{aligned} 
\qquad (2) $$

where with ⊗ being the tensor product,

$$ A = (\mathbf{I}_L \otimes \mathcal{F}_{\Omega}) \mathcal{S}, \quad \mathcal{S} = \begin{bmatrix} \mathbf{S}_1 \\ \vdots \\ \mathbf{S}_L \end{bmatrix}, \quad \mathbf{y} = \begin{bmatrix} \mathbf{y}_1 \\ \vdots \\ \mathbf{y}_L \end{bmatrix}. $$

The NC-PDNet architecture is made up of multiple blocks or iterations as shown in 1. We carry a buffer $x_b$ through iterations rather than relying on a single estimate. In each block $i$, we first compute a data consistency (DC) term in the k-space as $y_{DC}^H = A x_b[0] - y$, followed by density compensation $x_{DC} = A^H dy_{DC}$. A refinement network is then applied in the image domain using a residual CNN with three consecutive convolutions and ReLU activation, represented by $x_{i+1} = x_i + f_{\theta_i}(x_b, x_{DC})$.

The model is trained end-to-end, with trainable layers of each block $f_{\theta_i}$ having distinct weights that are not shared across iterations. For this study, we used an optimized configuration for NC-PDNet to fit the model on a single GPU, with six unrolled iterations, a buffer size of 2, and 16 convolution filters.

<figure><img src="image_3.png" /><figcaption>Fig. 1. NC-PDNet architecture [7] for multi-coil MRI reconstruction, where $N_f$ is the number of convolutional filters, $N_P$ is the buffer size, and $N_C$ is the number of unrolled iterations.</figcaption></figure>

## 2.2. Data

In this work, we used raw data from the Calgary-Campinas dataset [8], which includes 167 3D T1-weighted gradient-recalled echo scans with 1 mm isotropic sagittal acquisitions collected on a clinical 3-T MRI scanner from healthy subjects. The data set contains 117 scans (70.0%) acquired with a 12-channel receiver coil and 50 scans (30.0%) with a 32-channel coil. The acquisition parameters were TR / TE / TI = 6.3ms / 2.6ms / 650ms (93 scans, 55.7%) or TR/TE/TI = 7.4ms / 3.1ms / 400ms (74 scans, 44.3%). The matrix size was $N_x \times N_y \times N_z = 256 \times 218 \times [170-180]$, with the slice-encoded direction (kz) partially sampled up to 85% of its matrix size.

Subsequently, retrospective under-sampled data were obtained by projecting multi-coil reference images onto non-Cartesian under-sampled trajectories using the NUFFT operator. The ground truth was the square root of the sum of squares of all multi-coil Cartesian-sampled images. We randomly split the data into training (70%), validation (10%), and test (20%) sets. To perform coil compression, we applied singular value decomposition to the raw k-space data, retaining the singular vectors that account for 99% of the cumulative explained variance. Accordingly, we reduced the 12-channel and 32-channel k-space data to four and seven compressed channels respectively.

## 2.3. Non-Cartesian Sampling Patterns

Different 3D trajectories were generated to match an acceleration factor of 6 (≈ $\frac{218 \times 170}{6}$ trajectories with 512 samples each), which results in scan time of nearly 1.5min for 1 mm isotropic whole brain MRI. For our comparison, a 3D radial trajectory with golden means-based structure [9] and a 3D cone trajectory [10], both with center-out shots, were used. Also, as the density of the radial trajectories drops drastically at higher frequencies, we compared the reconstructed image quality using twisted projection imaging (TPI) [11]. Finally, we also tested the reconstructions on the GoLF-SPARKLING [12] sampling trajectory, which extends the original SPARKLING approach [13] that was originally tested for T2*-weighted imaging. In GoLF-SPARKLING, each sampling trajectory passes through the center of k-space as a Cartesian line, resulting in Cartesian sampling at the center of k-space, shown in green in Fig. 2. In practice, the acceleration factor of the GoLF-SPARKLING trajectory was slightly higher (6.7) and had almost 42% of the center of the k space sampled in a Cartesian way.

<figure><img src="image_4.png" /><figcaption>Fig. 2. Trajectories are shown with a reduced number of shots for clarity. GoLF-SPARKLING combines a non-Cartesian SPARKLING portion (blue) with a grid-sampled low-frequency region (green).</figcaption></figure>

## 2.4. Implementation, training and evaluation details

All the code was implemented in PyTorch. Alternating between the k-space and image domains was achieved using forward and ad joint NUFFT operators. To address memory constraints, we used
the MRI-NUFFT package [14] as its NUFFT implementations are memory and computationally efficient and support PyTorch tensors for the CufiNufft and GpuNUFFT backends.

We used the mean absolute error (MAE) as the loss function, the Adam optimizer with a learning rate of $1e-3$, and a reduce-on-plateau scheduler. The training was carried out for 200 epochs for all the networks on a single A100 GPU with 80GB of VRAM. We evaluated the model's performance on the test set using the checkpoint giving the best validation scores.

To assess the impact of channel configuration, we first conducted two separate training sessions for each non-Cartesian undersampling pattern: one with only 12-channel data and the other with 32-channel data. We refer to this approach as *channel-specific* training. In addition, we conducted training on the entire dataset, which is designed as *channel-agnostic* training.

The metrics used for evaluating and comparing the networks are the peak signal-to-noise ratio (PSNR) and structural similarity index (SSIM).

# 3. RESULTS AND ANALYSIS

## 3.1. Channel-Specific vs. Channel-Agnostic Training

Box plots in Fig. 3 illustrate that channel-agnostic models, when evaluated on the 32-channel test subset, outperform the 32-channel-specific model trained and tested exclusively on 32-channel data, showing higher PSNR and SSIM scores. Meanwhile, the channel-agnostic model achieves comparable or slightly lower scores when tested against the 12-channel-specific model on 12-channel data. These findings suggest that training a single channel-agnostic model across all data is more practical and efficient, as it provides near-equivalent performance across both 32- and 12-channel configurations without the need for separate, channel-specific models.

<table><caption>**Table 1.** Mean PSNR/SSIM scores of NC-PDNet with channel-agnostic training on 12-channel and 32-channel test subsets.</caption><thead><tr><th>Trajectory</th><th>12-ch. PSNR/SSIM</th><th>32-ch. PSNR/SSIM</th></tr></thead><tbody><tr><td>TPI</td><td>33.03 / 0.932</td><td>40.03 / 0.972</td></tr><tr><td>3D Radial</td><td>33.65 / 0.931</td><td>40.49 / 0.976</td></tr><tr><td>3D Cones</td><td>33.82 / 0.934</td><td>40.55 / 0.977</td></tr><tr><td>GS w/CC</td><td>36.83 / 0.960</td><td>42.89 / 0.987</td></tr><tr><td>GS</td><td><strong>37.62 / 0.962</strong></td><td><strong>43.33 / 0.987</strong></td></tr></tbody></table>

## 3.2. Performance benchmarking of selected non-Cartesian sampling patterns

**Quantitative results.** Results in Tab. 1 show that the undersampling with the GOLF-SPARKLING (denoted GS in the table, with or without coil compression – w/CC) outperforms all other methods, while the 3D cone trajectory performs slightly better than the 3D radial and TPI patterns. GoLF-SPARKLING achieves higher PSNR and SSIM scores due to its design, which incorporates Cartesian sampling in the central region of k-space. This sampling pattern enables a more accurate estimation of sensitivity maps, extracted by applying the adjoint NUFFT to a density-compensated low-frequency central region of k-space. As a result, GoLF-SPARKLING allows us to achieve superior performance compared to fully non-Cartesian sampling patterns.

<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Quantitative results of NC-PDNet with different trajectories, ordered from left to right: TPI, 3D Radial, 3D Cones, GoLF-SPARKLING (GS w/CC), and GoLF-SPARKLING without coil compression (GS).</figcaption></figure>

**Qualitative results.** Visual evaluation of reconstructed MR images confirms the quantitative metrics. As illustrated in Fig. 4, the fine structures of the cerebellum are sharper and more accurately reconstructed using the GoLF-SPARKLING and TPI trajectories, in contrast to the blurriness observed with the radial and cone patterns. Additionally, GoLF-SPARKLING demonstrates a higher fidelity in recovering brain structures than TPI, which displays a darkened appearance in some areas relative to the ground truth. This underestimation of signal intensity in TPI could obscure finer details in specific regions, potentially affecting diagnostic interpretation.

## 3.3. Assessing the Impact of Coil Compression

Plots of SSIM and PSNR test scores over the cumulative training time (Fig. 5) show that, for the same training duration, NC-PDNet with GS channel-agnostic training without coil compression achieves an equivalent average SSIM score and an increase of approximately +0.45dB in average PSNR compared to the same training performed with coil compressed k-space input data. In addition, coil compression training is less stable, as illustrated by the significant drop in performance in Fig. 5. However, Fig. 4 shows a very close perceived visual quality between the reconstructed images of models trained with and without coil compression, suggesting that in case of memory constraints, coil compression is a key strategy to maintain good performances while allowing for scalability with respect to the dimension of the multi-coil input data set.

Enabling coil compression reduces memory footprint and computational demand as it is shown in Tab. 2, however training for more epochs is required to attain performance comparable to a model trained without coil compression. Importantly, coil compression does not compromise the structural integrity of the reconstructed images, although it may lead to a decrease in average PSNR. Additionally, using coil compression accelerates inference by approximately
<figure><img src="image_6.png" /><figcaption>**Fig. 4.** Reconstruction results of the 90th slice of file e14079s3_P09216.7 from the test set. The top row shows reconstructions by different methods, while the bottom row displays zoomed-in regions outlined by red frames. Volume-wise PSNR and SSIM scores are indicated in the top left corner of each image.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>(a) Mean Test PSNR Score vs. Cumulative Training Time</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>(b) Mean Test SSIM Score vs. Cumulative Training Time</figcaption></figure>

<caption>**Fig. 5.** Evolution of test scores over cumulative training time for channel-agnostic NC-PDNet training.</caption>

2.66 times and 4.26 times for 12-channel and 32-channel volumes, respectively (see Tab. 3).

## 4. DISCUSSION AND CONCLUSION

In this work, we extend NC-PDNet, a density-compensated unrolled neural network, to support the reconstruction of multi-coil 3D non-Cartesian k-space data. Our retrospective comparative study demonstrated that the unique design of the GoLF-SPARKLING sampling trajectories enables the recovery of higher-quality images

<table><caption>**Table 2.** Comparison of NC-PDNet with GS, trained with and without coil compression (w/CC and w/o CC, resp.), showing average epoch training time and peak GPU memory usage per training step. Training data: 116 volumes (69.8% 12-channel, 30.2% 32-channel).</caption><thead><tr><td>Data set</td><td>Training Time per Epoch (s)</td><td>Peak GPU Memory (GB)</td></tr></thead><tbody><tr><td>w/CC</td><td>1922</td><td>17.86</td></tr><tr><td>w/o CC</td><td>7416</td><td>47.54</td></tr></tbody></table>

<table><caption>**Table 3.** Comparison of Inference time(s) / GPU memory usage(GB) between NC-PDNet with and without coil compression (w/CC and w/o CC, resp.) for 12-ch. and 32-ch. test volumes.</caption><thead><tr><td></td><td>12-ch. volume</td><td>32-ch. volume</td></tr></thead><tbody><tr><td>w/CC</td><td>2.91(s) / 4.39(GB)</td><td>4.95(s) / 5.49(GB)</td></tr><tr><td>w/o CC</td><td>7.75(s) / 7.57(GB)</td><td>21.1(s) / 19.74(GB)</td></tr></tbody></table>

compared to other sampling patterns. Presenting a benchmark of GoLF-SPARKLING's performance against commonly used non-Cartesian readouts in the literature is essential before advancing to further experiments. Our results indicated that while disabling coil compression provides enhanced performances in PSNR, this improvement has little impact on final image quality. We conclude that coil compression training in a channel-agnostic set-up is both practical and memory-efficient, achieving performance levels comparable to more complex training configurations. We reconstruct a 6x accelerated 1 mm isotropic 32-channel whole-brain 3D image in 4.95 seconds with a GPU memory usage of just 5.49 GB, demonstrating the clinical viability of this approach. A key future direction is to validate NC-PDNet's performance on prospectively undersampled GoLF-SPARKLING data and assess its scalability in more challenging imaging setups, notably at higher resolution.
## 5. COMPLIANCE WITH ETHICAL STANDARDS

This study was conducted retrospectively using human subject data made available in open access [8]. Ethical approval was not required as confirmed by the license attached with the open access data.

## 6. ACKNOWLEDGMENTS

This work was granted access to the CCRT High-Performance Computing (HPC) facility under the Grant CCRT2024-tanaasma awarded by the Fundamental Research Division (DRF) of CEA. This work was granted access to IDRIS' HPC resources under the allocation 2023 AD011011153R3 made by GENCI. The concepts and information presented in this abstract are based on research results that are not commercially available. Future availability cannot be guaranteed.

## 7. REFERENCES

[1] Michael Lustig et al., "Sparse MRI: The application of compressed sensing for rapid MR imaging," vol. 58, no. 6, pp. 1182–1195, 2007.

[2] Katherine L Wright et al., "Non-cartesian parallel imaging reconstruction," *Journal of Magnetic Resonance Imaging*, vol. 40, no. 5, pp. 1022–1040, 2014.

[3] Kerstin Hammernik, Teresa Klatzer, Erich Kobler, Michael P Recht, Daniel K Sodickson, Thomas Pock, and Florian Knoll, "Learning a variational network for reconstruction of accelerated mri data," *Magnetic resonance in medicine*, vol. 79, no. 6, pp. 3055–3071, 2018.

[4] Gushan Zeng, Yi Guo, Jiaying Zhan, Zi Wang, Zongying Lai, Xiaofeng Du, Xiaobo Qu, and Di Guo, "A review on deep learning MRI reconstruction without fully sampled k-space,"

[5] Jonas Adler and Ozan Öktem, "Learned primal-dual reconstruction," *IEEE transactions on medical imaging*, vol. 37, no. 6, pp. 1322–1332, 2018.

[6] Reinhard Heckel, Mathews Jacob, Akshay Chaudhari, Or Perlman, and Efrat Shimron, "Deep learning for accelerated and robust MRI reconstruction," *Magnetic Resonance Materials in Physics, Biology and Medicine*, vol. 37, no. 3, pp. 335–368, 2024.

[7] Zaccharie Ramzi, GR Chaithya, Jean-Luc Starck, and Philippe Ciuciu, "NC-PDNet: A density-compensated unrolled network for 2D and 3D non-Cartesian MRI reconstruction," *IEEE Transactions on Medical Imaging*, vol. 41, no. 7, pp. 1625–1638, 2022.

[8] Youssef Beauferris, Jonas Teuwen, Dimitrios Karkalousos, Nikita Moriakov, Matthan Caan, George Yiasemis, Lívia Rodrigues, Alexandre Lopes, Helio Pedrini, Letícia Rittner, et al., "Multi-coil mri reconstruction challenge—assessing brain mri reconstruction models and their generalizability to varying coil configurations," *Frontiers in neuroscience*, vol. 16, pp. 919186, 2022.

[9] Rachel W Chan, Elizabeth A Ramsay, Charles H Cunningham, and Donald B Plewes, "Temporal stability of adaptive 3d radial mri using multidimensional golden means," *Magnetic Resonance in Medicine: An Official Journal of the International Society for Magnetic Resonance in Medicine*, vol. 61, no. 2, pp. 354–363, 2009.

[10] BW Clare and DL Kepert, "The optimal packing of circles on a sphere," *Journal of mathematical chemistry*, vol. 6, no. 1, pp. 325–349, 1991.

[11] Fernando E. Boada, Joseph S. Gillen, Gary X. Shen, Sam Y Chang, and Keith R. Thulborn, "Fast three dimensional sodium imaging," vol. 37, no. 5, pp. 706–715, 1997.

[12] Chaithya Giliyar Radhakrishna, Guillaume Daval-Frérot, Aurélien Massire, Alexandre Vignaud, and Philippe Ciuciu, "Improving spreading projection algorithm for rapid k-space sampling trajectories through minimized off-resonance effects and gridding of low frequencies," *Magnetic Resonance in Medicine*, vol. 90, no. 3, pp. 1069–1085, 2023.

[13] G. R. Chaithya, Pierre Weiss, Guillaume Daval-Frérot, Aurélien Massire, Alexandre Vignaud, and Philippe Ciuciu, "Optimizing Full 3D SPARKLING Trajectories for High-Resolution Magnetic Resonance Imaging," *IEEE Transactions on Medical Imaging*, vol. 41, no. 8, pp. 2105–2117.

[14] P.-A. Comby, Guillaume Daval-Frerot, Chaithya Gr, A. Vignaud, and P. Ciuciu, "MRI-NUFFT: An open source Python package to make non-Cartesian MR Imaging easier," in *ISMRM Annual Meeting, (in Press)*.