Online MR Image Reconstruction for Compressed Sensing Acquisition in T2* Imaging P. Ciuciu1, 2 Joint work with: L. El Gueddari1,2, E. Chouzenoux3 and J-C. Pesquet3 1CEA/NeuroSpin, Univ. Paris-Saclay, Gif-sur-Yvette, France 2Inria, Parietal team, Univ. Paris-Saclay, France 3Inria, Opis team, CVN, Centrale-Supélec, Univ. Paris-Saclay, France October, 9 2020 L2S – CNRS – Supélec – Université Paris-Saclay, France P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 1 / 50

Outline 1 Context & Motivation Ultra High Field Imaging High resolution MRI Sampling in MRI 2 Multi-channel image reconstruction for non-Cartesian acquisition State-of the art Calibrationless reconstruction: Playing with the regularization Experiments & Results Summary 3 Online CS-MRI reconstruction Problem statement: The single channel case Multi-channel online reconstruction Experiments & Results Summary 4 Conclusion & Outlook P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 2 / 50

Outline 1 Context & Motivation Ultra High Field Imaging High resolution MRI Sampling in MRI 2 Multi-channel image reconstruction for non-Cartesian acquisition State-of the art Calibrationless reconstruction: Playing with the regularization Experiments & Results Summary 3 Online CS-MRI reconstruction Problem statement: The single channel case Multi-channel online reconstruction Experiments & Results Summary 4 Conclusion & Outlook P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 3 / 50

High resolution imaging: Increasing the SNR Ultra-high magnetic eld Signal-to-Noise Ratio (SNR) increases with the eld strength: SNR ∝ B1.65 0 1 Particularly suited for the T∗ 2 -weighted imaging contrast 1Pohmann, Speck, and Sche er 2016, Magnetic resonance in medicine. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 4 / 50

High resolution imaging: Increasing the SNR Ultra-high magnetic eld Signal-to-Noise Ratio (SNR) increases with the eld strength: SNR ∝ B1.65 0 Particularly suited for the T∗ 2 -weighted imaging contrast1 2017: 7 Tesla MR systems FDA-cleared & CE-marked for clinical practice 1Zwanenburg et al. 2010, European radiology. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 4 / 50

High resolution imaging: Increasing the SNR Multi-channel array coils Parallel imaging: collect multiple k-space data using a multi-receiver coil to boost the SNR 2. Illustration of multi-receiver coil (phased array). 2Roemer et al. 1990, Magnetic Resonance in Medicine P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 5 / 50

High resolution imaging: Increasing the SNR Multi-channel array coils Parallel imaging: collect multiple k-space data using a multi-receiver coil to boost the SNR 2. Ex-vivo baboon brain acquired at 7T with a single channel and 32-channel receiver coils. Single channel acquisition. Multi-channel acquisition. Courtesy Carole Lazarus 2Roemer et al. 1990, Magnetic Resonance in Medicine P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 5 / 50

Highest in vivo Human Brain data acquired at 7T Towards imaging at the meso-scale ( 100µm2 in plane) Highest in-vivo Human Brain MRI data collected at 7T 3: • B0 eld strength: 7T • Multi-channel receiver coil: 32 Rx • Resolution: 0.12 × 0.12 × 0.6mm3 • Field of View: 20.28 × 20.93 × 1.26 cm3 • Averages: 2 • Total scan time: 50 min 3Stucht et al. 2015, PloS one P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 6 / 50

Highest in vivo Human Brain data acquired at 7T Towards imaging at the meso-scale ( 100µm2 in plane) Highest in-vivo Human Brain MRI data collected at 7T 3: • B0 eld strength: 7T • Multi-channel receiver coil: 32 Rx • Resolution: 0.12 × 0.12 × 0.6mm3 • Field of View: 20.28 × 20.93 × 1.26 cm3 • Averages: 2 • Total scan time: 50 min How can we speed-up the acquisition? 3Stucht et al. 2015, PloS one P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 6 / 50

Acquisition in MRI Data collection in MRI is typically performed over a Cartesian grid: P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 7 / 50

Acquisition in MRI Data collection in MRI is typically performed over a Cartesian grid: P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 7 / 50

Acquisition in MRI Data collection in MRI is typically performed over a Cartesian grid: P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 7 / 50

Acquisition in MRI Data collection in MRI is typically performed over a Cartesian grid: P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 7 / 50

Compressed Sensing in MRI: Down-sample the k-space Speeding up the acquisition P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 8 / 50

Compressed Sensing: Beyond Nyquist sampling rate Sparsity Compressed sensing4,5 (CS) relies on three main assumptions: 1. Sparsity: Image to be reconstructed is represented by only a few non-zero coe cients in an transformed domain. 4Candès, Romberg, and Tao 2006, IEEE Transactions on information theory. 5Lustig, Donoho, and Pauly 2007, Magnetic Resonance in Medicine. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 9 / 50

Compressed Sensing: Beyond Nyquist sampling rate Incoherence Compressed sensing6,7 (CS) relies on three main assumptions: 1. Sparsity 2. Incoherence (local coherence8: µ(A) = n k=1 ak 2 ∞ , with A = F Ψ) Exact recovery with probability 1 − ε if m ≥ Cµ(A)s log( n ε ) 6Candès, Romberg, and Tao 2006, IEEE Transactions on information theory. 7Lustig, Donoho, and Pauly 2007, Magnetic Resonance in Medicine. 8Puy, Vandergheynst, and Wiaux 2011, IEEE Signal Processing Letters; Chau fert et al. 2014, SIAM Journal on Imaging Sciences. 8Lazarus et al. 2019, Magnetic Resonance in Medicine P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 10 / 50

Compressed Sensing: Beyond Nyquist sampling rate Incoherence Compressed sensing6,7 (CS) relies on three main assumptions: 1. Sparsity 2. Incoherence (local coherence8: µ(A) = n k=1 ak 2 ∞ , with A = F Ψ) Exact recovery with probability 1 − ε if m ≥ Cµ(A)s log( n ε ) Sparkling VDS9 Radial VDS 6Candès, Romberg, and Tao 2006, IEEE Transactions on information theory. 7Lustig, Donoho, and Pauly 2007, Magnetic Resonance in Medicine. 8Puy, Vandergheynst, and Wiaux 2011, IEEE Signal Processing Letters; Chau fert et al. 2014, SIAM Journal on Imaging Sciences. 9Lazarus et al. 2019, Magnetic Resonance in Medicine P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 10 / 50

Compressed Sensing: Beyond Nyquist sampling rate Incoherence Compressed sensing6,7 (CS) relies on three main assumptions: 1. Sparsity 2. Incoherence (local coherence8: µ(A) = n k=1 ak 2 ∞ , with A = F Ψ) Exact recovery with probability 1 − ε if m ≥ Cµ(A)s log( n ε ) Sparkling VDS9 Radial VDS 6Candès, Romberg, and Tao 2006, IEEE Transactions on information theory. 7Lustig, Donoho, and Pauly 2007, Magnetic Resonance in Medicine. 8Puy, Vandergheynst, and Wiaux 2011, IEEE Signal Processing Letters; Chau fert et al. 2014, SIAM Journal on Imaging Sciences. 9Lazarus et al. 2019, Magnetic Resonance in Medicine P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 10 / 50

Compressed Sensing: Beyond Nyquist sampling rate Nonlinear reconstruction: when sparsity meets incoherence 1. Sparsity 2. Local incoherence 3. Nonlinear sparsity-promoting reconstruction: x = arg min x∈CN 1 2 FΩS x − yΩS 2 F + λ Ψx 1, (1) with: • S: total number of shots acquired (with C samples each) • ΩS: support of under-sampling scheme in k-space • FΩS : under-sampled Fourier operator (e.g. NUFFTa for non-Cartesian sampling) • yΩS : complex-values k-space measurements (S × C N) • Ψ: sparsifying transform (e.g. wavelet transform (WT)) • λ > 0: Regularization parameter • x: MR image estimate P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 11 / 50

Compressed Sensing: O ine reconstruction First scan then reconstruct Two-step scan procedure: First: Collect data Second: Perform image recon P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 12 / 50

Compressed Sensing: Beyond Nyquist sampling rate Successful application of Compressed Sensing in MRI10 Sparkling VDS Radial VDS 9Lazarus et al. 2019, Magnetic Resonance in Medicine P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 13 / 50

Compressed Sensing: Beyond Nyquist sampling rate Successful application of Compressed Sensing in MRI10 Sparkling VDS Radial VDS Reconstruction takes 1-2 min for a 2D slice and the 3D scan consists of 208 slices! 9Lazarus et al. 2019, Magnetic Resonance in Medicine P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 13 / 50

Outline 1 Context & Motivation Ultra High Field Imaging High resolution MRI Sampling in MRI 2 Multi-channel image reconstruction for non-Cartesian acquisition State-of the art Calibrationless reconstruction: Playing with the regularization Experiments & Results Summary 3 Online CS-MRI reconstruction Problem statement: The single channel case Multi-channel online reconstruction Experiments & Results Summary 4 Conclusion & Outlook P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 14 / 50

Non-Cartesian multi-channel MR image reconstruction State of the art method Methods can be split in two categories: 1. Recovering a full FOV image: Self-Calibrated reconstruction 2. Recovering an image per channel: Calibrationless reconstruction P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 15 / 50

Non-Cartesian multi-channel MR image reconstruction Self-Calibrated reconstruction 1. Self-Calibrated reconstruction: • Solve an inverse problem and recover a single full FOV image: x = arg min x∈CN 1 2 L =1 σ−2 FΩ S x − y 2 2 + λ Ψx 1 (2) • Model the coil sensitivity pro le S for all channels = 1, . . . , L11,12 11Samsonov et al. 2004, Magnetic Resonance in Medicine. 12Uecker et al. 2014, Magnetic Resonance in Medicine. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 16 / 50

Non-Cartesian multi-channel MR image reconstruction Self-Calibrated reconstruction 1. Self-Calibrated reconstruction: • Solve an inverse problem and recover a single full FOV image: x = arg min x∈CN 1 2 L =1 σ−2 FΩ S x − y 2 2 + λ Ψx 1 (2) • Model the coil sensitivity pro le S for all channels = 1, . . . , L11,12 • Coil sensitivity pro les are subject/scan-speci c → estimated/extracted for each scan 11Samsonov et al. 2004, Magnetic Resonance in Medicine. 12Uecker et al. 2014, Magnetic Resonance in Medicine. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 16 / 50

Non-Cartesian multi-channel MR image reconstruction Calibrationless reconstruction 1. Self-calibrated methods, pros & cons: • good reconstruction performances • sensitivity pro les need to be set beforehand =⇒ pre-scan lengthen the acquisition • estimation must be performed during the acquisition =⇒ lengthen the reconstruction P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 17 / 50

Non-Cartesian multi-channel MR image reconstruction Calibrationless reconstruction 1. Self-calibrated methods, pros & cons: • good reconstruction performances • sensitivity pro les need to be set beforehand =⇒ pre-scan lengthen the acquisition • estimation must be performed during the acquisition =⇒ lengthen the reconstruction 2. Calibrationless reconstruction • No longer need coil sensitivity pro les • Exploit redundancy across channels to impose structured sparsity P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 17 / 50

Problem statement Calibration-less MR image reconstruction problem solved using an analysis formulation: X = arg min X∈CN×L 1 2 FΩ X − Y 2 2 + R(ΨX). (3) with: • Y = [y1, · · · , yL ] with y ∈ CM the th channel-speci c k-space • X = [x1, · · · , xL ] with x ∈ CN the th channel-speci c reconstructed image. • FΩ is the forward under-sampling Fourier operator • Ψ ∈ CNΨ×N linear operator related to a sparse decomposition • R is a convex regularization term that promotes sparsity with an explicit proximity operator. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 18 / 50

Optimization algorithm Primal dual optimization We aim to nd: X ∈ argmin X∈CN×L [f (X) + R(ΨX)] (4) where: f is convex, di ferentiable on CN×L and its gradient is β-Lipschitz R ∈ Γ0 (CNΨ×L) with a closed form proximity operator13, given by: prox R (Z) = argmin V ∈CNΨ×L 1 2 Z − V 2 + R(V ) (5) Note: Those are standard assumptions in optimization-based image reconstruction methods. 13Moreau 1962, Comptes Rendus de l’Académie des Sciences de Paris. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 19 / 50

Optimization algorithm Condat-V˜ u sequence Using a primal-dual optimization method proposed by Condat14-V˜ u15: Algorithm 1: Condat-V˜ u algorithm Data: X0, Z0 Result: XK , ZK initialize k = 0, τ > 0, κ > 0,; while k ≤ K do Xk+1 := Xk − τ ∇f (Xk ) + Ψ∗Zk ; Wk+1 := Zk + κ Ψ (2Xk+1 − Xk ); Zk+1 := Wk+1 − κ proxg/κ Wk+1 κ ; • if 1 τ − κ|||Ψ|||2 ≥ β 2 then the algorithm weakly converges to the solution of Eq. (4). • τ and κ hyper-parameters set as follows: τ := 1 β , κ := β 2|||Ψ|||2 14Condat 2013, Journal of Optimization Theory and Applications. 15V˜ u 2013, Advances in computational mathematics. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 20 / 50

Calibrationless reconstruction: playing with clustering regularization Clustering using the 2-norm Clustering using the ∞-norm P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 21 / 50

Regularizing through the Group-LASSO Group-LASSO (GL)16 regularization has been used in multi-task learning in di ferent eld including MRI reconstruction17. Ψ = 16Yuan and Lin 2006, Journal of the Royal Statistical Society: Series B (Statistical Methodology). 17Majumdar and Ward 2012, Magnetic Resonance in Medicine. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 22 / 50

Regularizing through the Group-LASSO Group-LASSO (GL)16 regularization has been used in multi-task learning in di ferent eld including MRI reconstruction17. The group-LASSO penalty is de ned as follows: RGL (Z) = g∈G Zg 2 • G is the set of groups de ning a partition of Z Z = 16Yuan and Lin 2006, Journal of the Royal Statistical Society: Series B (Statistical Methodology). 17Majumdar and Ward 2012, Magnetic Resonance in Medicine. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 22 / 50

Regularizing through the Group-LASSO Group-LASSO (GL)16 regularization has been used in multi-task learning in di ferent eld including MRI reconstruction17. The group-LASSO penalty is de ned as follows: RGL (Z) = Z 2,1 = λ NΨ p=1 L =1 |zp |2 with: • Z = ΨX • λ is a positive hyper-parameter Z = 16Yuan and Lin 2006, Journal of the Royal Statistical Society: Series B (Statistical Methodology). 17Majumdar and Ward 2012, Magnetic Resonance in Medicine. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 22 / 50

Regularizing through the Group-LASSO Group-LASSO (GL)16 regularization has been used in multi-task learning in di ferent eld including MRI reconstruction17. The group-LASSO penalty is de ned as follows: RGL (Z) = Z 2,1 = λ NΨ p=1 L =1 |zp |2 Group-LASSO regularization provides tighter recovery guarantees18 Z = 16Yuan and Lin 2006, Journal of the Royal Statistical Society: Series B (Statistical Methodology). 17Majumdar and Ward 2012, Magnetic Resonance in Medicine. 18Chun, Adcock, and Talavage 2016, IEEE Transactions on Medical Imaging P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 22 / 50

Joint sparsity regularization Octagonal Shrinkage and Clustering Algorithm for Regression Inferring the structure via a pairwise ∞ norm. OSCAR regularization19: Rc-OSCAR (Z) = NΨ p=1 λ L =1 |zpj | + γ

Experimental set-up Sequence parameters: • Ex-vivo human brain • 7T Siemens Scanner • 1Tx/32Rx Nova coil • GRE Sparkling trajectory • 0.390 x 0.390 x 1.5 m3 resolution • Acceleration factor of 20 in time • Under-sampling factor of 2.5 • Ψ: Daubechies 4 transform Figure: Sparkling trajectory Figure: Ex-vivo Human brain Hyper-parameters set using a grid-search procedure to maximize the SSIM score. Cartesian scan 512×512 was acquired and used for reference. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 24 / 50

Results Similarity score iFT 1-ESPIRIT20 p-LORAKS21 group-LASSO22 c-OSCAR23 k-support norm SSIM 0.884 0.885 0.753 0.897 0.901 0.900 pSNR 28.25 26.48 25.52 28.59 29.77 30.29 NRMSE 0.1911 0.2276 0.2536 0.1859 0.1604 0.1510 20Uecker et al. 2014, Magnetic Resonance in Medicine. 21Haldar and Zhuo 2016, Magnetic resonance in medicine. 22Majumdar and Ward 2012, Magnetic Resonance in Medicine. 23El Gueddari et al. 2019, ISBI 2019. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 25 / 50

Results Magnitude images Reference 1-ESPIRIT p-LORAKS group-LASSO c-OSCAR k-support norm P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 26 / 50

Results Magnitude images Reference 1-ESPIRIT p-LORAKS group-LASSO c-OSCAR k-support norm P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 26 / 50

Interim summary Summary: • Investigate two classes of MR image reconstruction methods for non-Cartesian parallel imaging acquisitions • Impose structured sparsity for calibrationless MR image reconstruction • Achieve same image quality in the calibrationless setting as in the self-calibrated one at the cost of a larger memory footprint and computation cost → How can we speed-up the reconstruction. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 27 / 50

Outline 1 Context & Motivation Ultra High Field Imaging High resolution MRI Sampling in MRI 2 Multi-channel image reconstruction for non-Cartesian acquisition State-of the art Calibrationless reconstruction: Playing with the regularization Experiments & Results Summary 3 Online CS-MRI reconstruction Problem statement: The single channel case Multi-channel online reconstruction Experiments & Results Summary 4 Conclusion & Outlook P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 28 / 50

Speeding-up the reconstruction Structural MRI acquisition is segmented in time P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 29 / 50

Speeding-up the reconstruction Structural MRI acquisition is segmented in time P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 29 / 50

Online CS-MRI reconstruction Single-channel receiver case Multi-shot acquisition: ∀i ∈ {1, . . . , S}, collect data yi over Γi support with |Γi | = C P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 30 / 50

Online CS-MRI reconstruction Single-channel receiver case Multi-shot acquisition: ∀i ∈ {1, . . . , S}, collect data yi over Γi support with |Γi | = C Denote ωk = {Γ(k−1)∗bs +1 , . . . , Γk∗bs } the bs-dimensional set of shots to be added to the kth mini-batch de ned as follows: Ωk = ∪k i=1 ωi , with k ≤ S/bs Solve partial image reconstruction: ∀k ∈ {1, . . . , S}, xk = arg min x∈CN S 2k FΩk x − yΩk 2 F + λ Ψx 1, (6) with: • FΩk : undersampled Fourier operator over Ωk • yΩk : complex-valued k-space measurements over Ωk • xk : estimated MR image from incomplete data yΩk Online time constraints: nb × Tit ≈ bs × TR. The problem is convex → nal solution does not depend on the initialization. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 30 / 50

Online optimization algorithm Single-channel receiver case To solve Pb (6), we adapted the Condat-V˜ u algorithm to the online framework. Algorithm 1: Online mini-batch reconstruction algorithm for solving Problem (6). 1 initialize k = bs , x0, v0; 2 while k ≤ S − bs do 3 τk = 1/βk ; 4 κk = βk /(2|||Ψ|||2); 5 for t = 1, 2, . . . , nb do 6 xk t = xk t−1 − τk ∇fΩk (xk t−1 ) + Ψ∗vk t−1 ; 7 wk t = vk t−1 + κk Ψ 2xk t − xk t−1 ; 8 vk t = wk t − κk proxg/κk wk t /κk ; 9 (xk+bs 0 , vk+bs 0 ) ← (xk nb , vk nb ); 10 k ← k + bs ; P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 31 / 50

Online optimization algorithm Single-channel receiver case To solve Pb (6), we adapted the Condat-V˜ u algorithm to the online framework. Algorithm 1: Online mini-batch reconstruction algorithm for solving Problem (6). 1 initialize k = bs , x0, v0; 2 while k ≤ S − bs do 3 τk = 1/βk ; 4 κk = βk /(2|||Ψ|||2); 5 for t = 1, 2, . . . , nb do 6 xk t = xk t−1 − τk ∇fΩk (xk t−1 ) + Ψ∗vk t−1 ; 7 wk t = vk t−1 + κk Ψ 2xk t − xk t−1 ; 8 vk t = wk t − κk proxg/κk wk t /κk ; 9 (xk+bs 0 , vk+bs 0 ) ← (xk nb , vk nb ); 10 k ← k + bs ; Primal variable P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 31 / 50

Online optimization algorithm Single-channel receiver case To solve Pb (6), we adapted the Condat-V˜ u algorithm to the online framework. Algorithm 1: Online mini-batch reconstruction algorithm for solving Problem (6). 1 initialize k = bs , x0, v0; 2 while k ≤ S − bs do 3 τk = 1/βk ; 4 κk = βk /(2|||Ψ|||2); 5 for t = 1, 2, . . . , nb do 6 xk t = xk t−1 − τk ∇fΩk (xk t−1 ) + Ψ∗vk t−1 ; 7 wk t = vk t−1 + κk Ψ 2xk t − xk t−1 ; 8 vk t = wk t − κk proxg/κk wk t /κk ; 9 (xk+bs 0 , vk+bs 0 ) ← (xk nb , vk nb ); 10 k ← k + bs ; Dual variable P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 31 / 50

Online optimization algorithm Single-channel receiver case To solve Pb (6), we adapted the Condat-V˜ u algorithm to the online framework. Algorithm 1: Online mini-batch reconstruction algorithm for solving Problem (6). 1 initialize k = bs , x0, v0; 2 while k ≤ S − bs do 3 τk = 1/βk ; 4 κk = βk /(2|||Ψ|||2); 5 for t = 1, 2, . . . , nb do 6 xk t = xk t−1 − τk ∇fΩk (xk t−1 ) + Ψ∗vk t−1 ; 7 wk t = vk t−1 + κk Ψ 2xk t − xk t−1 ; 8 vk t = wk t − κk proxg/κk wk t /κk ; 9 (xk+bs 0 , vk+bs 0 ) ← (xk nb , vk nb ); 10 k ← k + bs ; Primal gradient step P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 31 / 50

Online optimization algorithm Single-channel receiver case To solve Pb (6), we adapted the Condat-V˜ u algorithm to the online framework. Algorithm 1: Online mini-batch reconstruction algorithm for solving Problem (6). 1 initialize k = bs , x0, v0; 2 while k ≤ S − bs do 3 τk = 1/βk ; 4 κk = βk /(2|||Ψ|||2); 5 for t = 1, 2, . . . , nb do 6 xk t = xk t−1 − τk ∇fΩk (xk t−1 ) + Ψ∗vk t−1 ; 7 wk t = vk t−1 + κk Ψ 2xk t − xk t−1 ; 8 vk t = wk t − κk proxg/κk wk t /κk ; 9 (xk+bs 0 , vk+bs 0 ) ← (xk nb , vk nb ); 10 k ← k + bs ; Dual proximal step P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 31 / 50

Online optimization algorithm Single-channel receiver case To solve Pb (6), we adapted the Condat-V˜ u algorithm to the online framework. Algorithm 1: Online mini-batch reconstruction algorithm for solving Problem (6). 1 initialize k = bs , x0, v0; 2 while k ≤ S − bs do 3 τk = 1/βk ; 4 κk = βk /(2|||Ψ|||2); 5 for t = 1, 2, . . . , nb do 6 xk t = xk t−1 − τk ∇fΩk (xk t−1 ) + Ψ∗vk t−1 ; 7 wk t = vk t−1 + κk Ψ 2xk t − xk t−1 ; 8 vk t = wk t − κk proxg/κk wk t /κk ; 9 (xk+bs 0 , vk+bs 0 ) ← (xk nb , vk nb ); 10 k ← k + bs ; Condat-V˜ u sequence P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 31 / 50

Online optimization algorithm Single-channel receiver case To solve Pb (6), we adapted the Condat-V˜ u algorithm to the online framework. Algorithm 1: Online mini-batch reconstruction algorithm for solving Problem (6). 1 initialize k = bs , x0, v0; 2 while k ≤ S − bs do 3 τk = 1/βk ; 4 κk = βk /(2|||Ψ|||2); 5 for t = 1, 2, . . . , nb do 6 xk t = xk t−1 − τk ∇fΩk (xk t−1 ) + Ψ∗vk t−1 ; 7 wk t = vk t−1 + κk Ψ 2xk t − xk t−1 ; 8 vk t = wk t − κk proxg/κk wk t /κk ; 9 (xk+bs 0 , vk+bs 0 ) ← (xk nb , vk nb ); 10 k ← k + bs ; Warm restart P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 31 / 50

Online optimization algorithm Single-channel receiver case To solve Pb (6), we adapted the Condat-V˜ u algorithm to the online framework. Algorithm 1: Online mini-batch reconstruction algorithm for solving Problem (6). 1 initialize k = bs , x0, v0; 2 while k ≤ S − bs do 3 τk = 1/βk ; 4 κk = βk /(2|||Ψ|||2); 5 for t = 1, 2, . . . , nb do 6 xk t = xk t−1 − τk ∇fΩk (xk t−1 ) + Ψ∗vk t−1 ; 7 wk t = vk t−1 + κk Ψ 2xk t − xk t−1 ; 8 vk t = wk t − κk proxg/κk wk t /κk ; 9 (xk+bs 0 , vk+bs 0 ) ← (xk nb , vk nb ); 10 k ← k + bs ; Concatenate batches P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 31 / 50

Acquisition parameters for prospective CS Single-channel receiver case Acquisition parameters • Field strength: 7 Tesla • Repetition Time (TR): 550 ms • Echo Time (TE): 30 ms • Field of View: 20.4 cm • Resolution: 0.4 × 0.4 × 3 mm3 • Matrix size: 512 × 512 • Trajectory: Sparkling26 • Number of shots (S): 34 • Acceleration/undersampling factors: 15/2.5 Cartesian Reference Figure: Ex-vivo baboon brain Non-Cartesian sampling Figure: 15-fold accelerated Sparkling trajectories 26Lazarus et al. 2019, Magnetic Resonance in Medicine P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 32 / 50

Reconstruction parameters for online image processing Single-channel receiver case Online processing constraints nb × Tit ≈ bs × TR. • Fourier Transform: GPU NUFFT27 • Time per iteration: Tit 78.2 ± 8.9 ms • Nb. shots/mini-batch: bs ∈ {1, 2, 17, 34} • Nb. iterations/mini-batch: nb ∈ {5, 11, 93, 200} Regularization term • Sparse transform: 2D decimated Symmlet-8 WT • Hyper-parameter: λ retrospectively tuned 27Lin 2018, Journal of Imaging P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 33 / 50

Non-Cartesian online CS-MR image reconstruction Single-channel receiver case Figure: Evolution of the global cost function (left) and SSIM score (right) over time. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 34 / 50

Non-Cartesian online CS-MR image reconstruction Single-channel receiver case Figure: Evolution of the global cost function (left) and SSIM score (right) over time. =⇒ The smaller the mini-batch size, the sooner and better the image solution. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 34 / 50

Non-Cartesian online CS-MR image reconstruction Single-channel receiver case Evolution of the output xk over the acquisition, mini-batch size of 1. Figure: Collected shots over time Figure: Evolution of xk over time P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 35 / 50

Online reconstruction Extension to the multi-channel acquisition Reconstruction method for multi-channel non-Cartesian acquisitions: 1. Self-Calibrated reconstruction Provides good reconstruction Requires the coil sensitivity maps estimation Gradient-Lipschitz constant depends on the sensitivity pro le estimation 2. Calibrationless reconstruction Slightly better reconstruction Does not rely on coil sensitivity pro les Lipschitz constant only depends on the sampling pattern → More suited for online reconstruction 27El Gueddari et al. 2020 P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 36 / 50

Online reconstruction Extension to the multi-channel acquisition Rely on calibrationless reconstruction: Xk = arg min X∈CN×L 1 2 FΩk X − YΩk 2 2 + R(ΨX). (7) with sub-band version of OSCAR regularization28: Z = ΨX, Rb-OSCAR (Z) = S s=1 B b=1 PbL =1 λ|zs,b, | + γ S s=1 B b=1 1≤ < ≤PbL max{|zs,b, |, |zs,b, |}, (8) λ > 0, γ > 0. 28El Gueddari et al. 2020 P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 36 / 50

Acquisition parameters for prospective CS Multi-channel receiver case Acquisition parameters • Field strength: 7 Tesla • Nb of coils (L): 32 • Repetition Time (TR): 550 ms • Field of View: 20.4 cm • Resolution: 0.4 × 0.4 × 1.5 mm3 • Trajectory: Sparkling29 • Number of spokes (S): 34 • Acceleration/undersampling factors: 15/2.5 Cartesian reference Figure: Ex-vivo human brain Non-Cartesian sampling Figure: 15-fold accelerated Sparkling trajectories 29Lazarus et al. 2019, Magnetic Resonance in Medicine P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 37 / 50

Reconstruction parameters for online image processing Multi-channel receiver case Online processing constraints nb × Tit ≈ bs × TR. • Fourier Transform: GPU NUFFT30 • Time per iteration: Tit 4.29 ± 0.111 s • Nb. shots/mini-batch: bs ∈ {17, 34} • Nb. iterations/mini-batch: nb ∈ {1, 200} Regularization term • Sparse transform: 2D decimated Symmlet-8 WT • Hyper-parameter: (λ, γ) retrospectively tuned 30Lin 2018, Journal of Imaging P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 38 / 50

Results: Prospective non-Cartesian reconstruction Online multi-channel reconstruction Figure: Evolution of the global cost function (left) and estimated image obtained by the end of acquisition. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 39 / 50

Pro ling computation time Online multi-channel reconstruction Which term contributes the most to the global cost function? P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 40 / 50

Pro ling computation time Online multi-channel reconstruction Which term contributes the most to the global cost function? Gradient step Proximity op. step Linear Operator Total time per iteration direct adjoint S = 34 750 ms ± 32.2 ms 847 ms ± 17.9 ms 998 ms ± 15.8 ms 667 ms ± 16.2 ms 4.29 s ± 111 ms P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 40 / 50

Pro ling computation time Online multi-channel reconstruction Which term contributes the most to the global cost function? • The contribution of the data delity term is predominant during the rst iterations. • The gradient step only takes a quarter of computing time per iteration. =⇒ During acquisition: minimizing only the smooth term. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 40 / 50

The online trick for parallel imaging Multi-channel receiver case Algorithm 1: A fast method for online CS+PI calibration-less reconstruction initialize k = bs , X0, Z0; while k ≤ S − bs do for t = 1, 2, . . . , nb do Xk t = Xk t−1 − ∇fΩk (Xk t−1 ) /βk ; Xk+bs 0 ← Xk nb ; k ← k + bs ; ZS 0 ← ΨXS 0 ; τ ← 1/βS ; κ = βS /(2|||Ψ|||2); for t = 1, 2, . . . , 200 do XS t = XS t−1 − τ ∇fΩS (XS t−1 ) + Ψ∗ZS t−1 ; W S t = ZS t−1 + κΨ 2XS t − XS t−1 ; ZS t = W S t − κ proxg/κ W S t /κ ; Adaptive Gradient step Solving calibrationless problem via Condat-V˜ u algorithm P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 41 / 50

Online calibration-less CS-PI reconstruction Multi-channel receiver case Minimizing only the data- delity term allowed us to reduce the mini-batch size to bs = 2 shots. With the online trick P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 42 / 50

Online calibration-less CS-PI reconstruction Multi-channel receiver case Comparison with the fully online framework First approach With the online trick P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 42 / 50

Online calibration-less CS-PI reconstruction Multi-channel receiver case Recovered solution at convergence Reference Calibrationless reconstruction P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 42 / 50

Online reconstruction Partial summary • An online reconstruction approach for 2D MR image reconstruction • A Calibrationless version to make this online framework compatible with multi-channel acquisition and reconstruction • Small batches allow to achieve high image quality by the end of acquisition. • Delivery of good images by the end of acquisitions • The bigger the computational resources (GPU), the better the image quality! P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 43 / 50

Open source software: A fruit of the COSMIC project Python package for image reconstruction: Python Sparse data Analysis Package (PySAP) • Free open source Python package • pysap-mri: dedicated plugin for MRI reconstruction 31 • Unit-tested & Continuous integration • pysap-tutorials: Examples, Tutorials • Documentation 31El Gueddari et al. 2020, ISMRM workshop on Data Sampling P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 44 / 50

Outline 1 Context & Motivation Ultra High Field Imaging High resolution MRI Sampling in MRI 2 Multi-channel image reconstruction for non-Cartesian acquisition State-of the art Calibrationless reconstruction: Playing with the regularization Experiments & Results Summary 3 Online CS-MRI reconstruction Problem statement: The single channel case Multi-channel online reconstruction Experiments & Results Summary 4 Conclusion & Outlook P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 45 / 50

Conclusion & Outlook Conclusion: • Propose a provably convergent online reconstruction framework • Implement a mini-batch formulation that allows to interleave acquisition & recon and favor small batch sizes • Propose many algorithmic accelerations: GPU NUFFT, warm-restart, online trick, density compensation • Provide a feedback during the scan • Share the code and implement online recon on the scanner (Gadgetron) P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 46 / 50

Conclusion & Outlook Conclusion: • Propose a provably convergent online reconstruction framework • Implement a mini-batch formulation that allows to interleave acquisition & recon and favor small batch sizes • Propose many algorithmic accelerations: GPU NUFFT, warm-restart, online trick, density compensation • Provide a feedback during the scan • Share the code and implement online recon on the scanner (Gadgetron) Outlook: • Faster optimization algorithms (greedy-FISTA32, POGM33) for orthogonal WT • Further parallelization needed for 3D imaging 32Liang and Schönlieb 2019, Convergence 33Taylor, Hendrickx, and Glineur 2017, SIAM Journal on Optimization P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 46 / 50

Bibliography I Bondell, H.D. and B.J. Reich (2008). “Simultaneous regression shrinkage, variable selection, and supervised clustering of predictors with OSCAR”. In: Biometrics 64.1, pp. 115–123. Candès, Emmanuel J., J. Romberg, and T. Tao (2006). “Robust uncertainty principles: exact signal reconstruction from highly incomplete frequency information”. In: IEEE Transactions on information theory 52.2, pp. 489–509. Chau fert, Nicolas et al. (2014). “Variable density sampling with continuous trajectories”. In: SIAM Journal on Imaging Sciences 7.4, pp. 1962–1992. Chun, I.Y., B. Adcock, and T.M. Talavage (2016). “E cient compressed sensing SENSE pMRI reconstruction with joint sparsity promotion”. In: IEEE Transactions on Medical Imaging 35.1, pp. 354–368. Condat, Laurent (2013). “A primal–dual splitting method for convex optimization involving Lipschitzian, proximable and linear composite terms”. In: Journal of Optimization Theory and Applications 158.2, pp. 460–479. El Gueddari, L. et al. (Jan. 2020). Calibration-less parallel imaging compressed sensing reconstruction based on OSCAR regularization. IEEE Transactions on Medical Imaging. Saclay, France: CEA/NeuroSpin & INRIA Saclay & CVN Centrale-Supélec. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 47 / 50

Bibliography II El Gueddari, Loubna et al. (2019). “Calibrationless OSCAR-based image reconstruction in compressed sensing parallel MRI”. In: ISBI 2019. El Gueddari, Loubna et al. (2020). “PySAP-MRI: a Python Package for MR Image Reconstruction”. In: ISMRM workshop on Data Sampling. Fessler, Je frey A and Bradley P Sutton (2003). “Nonuniform fast Fourier transforms using min-max interpolation”. In: IEEE transactions on signal processing 51.2, pp. 560–574. Haldar, Justin P and Jingwei Zhuo (2016). “P-LORAKS: Low-rank modeling of local k-space neighborhoods with parallel imaging data”. In: Magnetic resonance in medicine 75.4, pp. 1499–1514. Lazarus, Carole et al. (2019). “SPARKLING: variable-density k-space lling curves for accelerated T2*-weighted MRI”. In: Magnetic Resonance in Medicine 81.6, pp. 3643–3661. Liang, Jingwei and Carola-Bibiane Schönlieb (2019). “Improving “Fast Iterative Shrinkage-Thresholding Algorithm”: Faster, Smarter and Greedier”. In: Convergence 50, p. 12. Lin, Jyh-Miin (2018). “Python Non-Uniform Fast Fourier Transform (PyNUFFT): An Accelerated Non-Cartesian MRI Package on a Heterogeneous Platform (CPU/GPU)”. In: Journal of Imaging 4.3, p. 51. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 48 / 50

Bibliography III Lustig, Miki, D. Donoho, and John M. Pauly (2007). “Sparse MRI: The application of compressed sensing for rapid MR imaging”. In: Magnetic Resonance in Medicine 58.6, pp. 1182–1195. Majumdar, Angshul and Rabab Ward (2012). “Calibration-less multi-coil MR image reconstruction”. In: Magnetic Resonance in Medicine 30.7, pp. 1032–1045. Moreau, Jean Jacques (1962). “Fonctions convexes duales et points proximaux dans un espace hilbertien”. In: Comptes Rendus de l’Académie des Sciences de Paris. Pohmann, Rolf, Oliver Speck, and Klaus Sche er (2016). “Signal-to-noise ratio and MR tissue parameters in human brain imaging at 3, 7, and 9.4 tesla using current receive coil arrays”. In: Magnetic resonance in medicine 75.2, pp. 801–809. Puy, Gilles, Pierre Vandergheynst, and Yves Wiaux (2011). “On variable density compressive sampling”. In: IEEE Signal Processing Letters 18.10, pp. 595–598. Roemer, P.B. et al. (1990). “The NMR phased array”. In: Magnetic Resonance in Medicine 16.2, pp. 192–225. Samsonov, Alexei A et al. (2004). “POCSENSE: POCS-based reconstruction for sensitivity encoded magnetic resonance imaging”. In: Magnetic Resonance in Medicine 52.6, pp. 1397–1406. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 49 / 50

Bibliography IV Stucht, Daniel et al. (2015). “Highest resolution in vivo human brain MRI using prospective motion correction”. In: PloS one 10.7, e0133921. Taylor, Adrien B, Julien M Hendrickx, and François Glineur (2017). “Exact worst-case performance of rst-order methods for composite convex optimization”. In: SIAM Journal on Optimization 27.3, pp. 1283–1313. Uecker, M. et al. (2014). “ESPIRiT– an eigenvalue approach to autocalibrating parallel MRI: where SENSE meets GRAPPA”. In: Magnetic Resonance in Medicine 71.3, pp. 990–1001. V˜ u, B.C (2013). “A splitting algorithm for dual monotone inclusions involving cocoercive operators”. In: Advances in computational mathematics 38.3, pp. 667–681. Yuan, Ming and Yi Lin (2006). “Model selection and estimation in regression with grouped variables”. In: Journal of the Royal Statistical Society: Series B (Statistical Methodology) 68.1, pp. 49–67. Zwanenburg, Jaco JM et al. (2010). “Fluid attenuated inversion recovery (FLAIR) MRI at 7.0 Tesla: comparison with 1.5 and 3.0 Tesla”. In: European radiology 20.4, pp. 915–922. P. Ciuciu (NeuroSpin) Online CS-MR image reconstruction 50 / 50