We propose a new approach for large-scale high-dynamic range computational imaging. Deep Neural Networks (DNNs) trained end-to-end can solve linear inverse imaging problems almost instantaneously. While unfolded architectures provide robustness to measurement setting variations, embedding large-scale measurement operators in DNN architectures is impractical. Alternative Plug-and-Play (PnP) approaches, where the denoising DNNs are blind to the measurement setting, have proven effective to address scalability and high-dynamic range challenges, but rely on highly iterative algorithms. We propose a residual DNN series approach, also interpretable as a learned version of matching pursuit, where the reconstructed image is a sum of residual images progressively increasing the dynamic range, and estimated iteratively by DNNs taking the back-projected data residual of the previous iteration as input. We demonstrate on radio-astronomical imaging simulations that a series of only few terms provides a reconstruction quality competitive with PnP, at a fraction of the cost.

Image restoration has long been one of the key research topics in image processing. Many mathematical approaches have been developed to solve this problem, e.g., variational methods, wavelet techniques, or Bayesian methods. With the widespread of neural network (NN) models in all the subdomains of data science, the performance limits of these methods are further pushed. One of the most successful strategies consists of plugging NNs in existing optimization algorithms. However, so doing raises several mathematical and practical challenges. One of the main issues is to secure the convergence of the resulting iterative scheme. Further questions concerning the characterization of the reached limit are also worth being addressed. In this paper, we show that the theory of maximally monotone operators allows us to bring insightful answers to these problems and to design firmly nonexpansive NNs; combining these with postprocessing NNs leads to excellent global restoration quality.

We introduce a new paradigm for solving regularized variational problems. These are typically formulated to address ill-posed inverse problems encountered in signal and image processing. The objective function is traditionally defined by adding a regularization function to a data fit term, which is subsequently minimized by using iterative optimization algorithms. Recently, several works have proposed to replace the operator related to the regularization by a more sophisticated denoiser. These approaches, known as plug-and-play (PnP) methods, have shown excellent performance. Although it has been noticed that, under some Lipschitz properties on the denoisers, the convergence of the resulting algorithm is guaranteed, little is known about characterizing the asymptotically delivered solution. In the current article, we propose to address this limitation. More specifically, instead of employing a functional regularization, we perform an operator regularization, where a maximally monotone operator (MMO) is learned in a supervised manner. This formulation is flexible as it allows the solution to be characterized through a broad range of variational inequalities, and it includes convex regularizations as special cases. From an algorithmic standpoint, the proposed approach consists in replacing the resolvent of the MMO by a neural network (NN). We present a universal approximation theorem proving that nonexpansive NNs are suitable models for the resolvent of a wide class of MMOs. The proposed approach thus provides a sound theoretical framework for analyzing the asymptotic behavior of first-order PnP algorithms. In addition, we propose a numerical strategy to train NNs corresponding to resolvents of MMOs. We apply our approach to image restoration problems and demonstrate its validity in terms of both convergence and quality.

During the last decades, reweighted procedures have shown high efficiency in computational imaging. They aim to handle non-convex composite penalization functions by iteratively solving multiple approximated sub-problems. Although the asymptotic behaviour of these methods has recently been investigated in several works, they all necessitate the sub-problems to be solved accurately, which can be sub-optimal in practice. In this work we present a reweighted forward-backward algorithm designed to handle non-convex composite functions. Unlike existing convergence studies in the literature, the weighting procedure is directly included within the iterations, avoiding the need for solving any sub-problem. We show that the obtained reweighted forward-backward algorithm converges to a critical point of the initial objective function. We illustrate the good behaviour of the proposed approach on a Fourier imaging example borrowed to radio-astronomical imaging.

Building nonexpansive Convolutional Neural Networks (CNNs) is a challenging problem that has recently gained a lot of attention from the image processing community. In particular, it appears to be the key to obtain convergent Plugand-Play algorithms. This problem, which relies on an accurate control of the the Lipschitz constant of the convolutional layers, has also been investigated for Generative Adversarial Networks to improve robustness to adversarial perturbations. However, to the best of our knowledge, no efficient method has been developed yet to build nonexpansive CNNs. In this paper, we develop an optimization algorithm that can be incorporated in the training of a network to ensure the nonexpansiveness of its convolutional layers. This is shown to allow us to build firmly nonexpansive CNNs. We apply the proposed approach to train a CNN for an image denoising task and show its effectiveness through simulations.

We propose a Bayesian uncertainty quantification method for large-scale imaging inverse problems. Our method applies to all Bayesian models that are log-concave, where maximum a posteriori (MAP) estimation is a convex optimization problem. The method is a framework to analyze the confidence in specific structures observed in MAP estimates (e.g., lesions in medical imaging, celestial sources in astronomical imaging), to enable using them as evidence to inform decisions and conclusions. Precisely, following Bayesian decision theory, we seek to assert the structures under scrutiny by performing a Bayesian hypothesis test that proceeds as follows: first, it postulates that the structures are not present in the true image, and then seeks to use the data and prior knowledge to reject this null hypothesis with high probability. Computing such tests for imaging problems is generally very difficult because of the high dimensionality involved. A main feature of this work is to leverage probability concentration phenomena and the underlying convex geometry to formulate the Bayesian hypothesis test as a convex problem, which we then efficiently solve by using scalable optimization algorithms. This allows scaling to high-resolution and high-sensitivity imaging problems that are computationally unaffordable for other Bayesian computation approaches. We illustrate our methodology, dubbed BUQO (Bayesian Uncertainty Quantification by Optimization), on a range of challenging Fourier imaging problems arising in astronomy and medicine. MATLAB code for the proposed uncertainty quantification method is available on GitHub.

We investigate the benefits of known partial support for the recovery of joint-sparse signals and demonstrate that it is advantageous in terms of recovery performance for both rank-blind and rank-aware algorithms. We suggest extensions of several joint-sparse recovery algorithms, e.g., simultaneous normalized iterative hard thresholding, subspace greedy methods and subspace-augmented multiple signal classification techniques. We describe a direct application of the proposed methods for compressive multiplexing of ultrasound (US) signals. The technique exploits the compressive multiplexer architecture for signal compression and relies on joint-sparsity of US signals in the frequency domain for signal reconstruction. We validate the proposed algorithms on numerical experiments and show their superiority against state-of-the-art approaches in rank-defective cases. We also demonstrate that the techniques lead to a significant increase of the image quality on in vivo carotid images compared to reconstruction without partially known support. The supporting code is available on https://github.com/AdriBesson/spl2018_joint_sparse.

This paper considers the problem of reconstructing raw signals from random projections in the context of time-of-flight imaging with an array of sensors. It presents a new signal model, coined as multi-channel pulse-stream model, which exploits pulse-stream models and accounts for additional structure induced by inter-sensor dependencies. We propose a sampling theorem and a reconstruction algorithm, based on ℓ -minimization, for signals belonging to such a model. We demonstrate the benefits of the proposed approach by means of numerical simulations and on a real non-destructive-evaluation application where the peak-signal-to-noise-ratio is increased by 3 dB compared to standard compressed-sensing strategies.

We develop techniques to solve ill-posed inverse problems on the sphere by sparse regularization, exploiting sparsity in both axisymmetric and directional scale-discretized wavelet space. Denoising, in painting, and deconvolution problems and combinations thereof, are considered as examples. Inverse problems are solved in both the analysis and synthesis settings, with a number of different sampling schemes. The most effective approach is that with the most restricted solution-space, which depends on the interplay between the adopted sampling scheme, the selection of the analysis/synthesis problem, and any weighting of the ℓ1norm appearing in the regularization problem. More efficient sampling schemes on the sphere improve reconstruction fidelity by restricting the solution-space and also by improving sparsity in wavelet space. We apply the technique to denoise Planck 353-GHz observations, improving the ability to extract the structure of Galactic dust emission, which is important for studying Galactic magnetism.

The theory of compressed sensing (CS) leverages upon structure of signals in order to reduce the number of samples needed to reconstruct a signal, compared to the Nyquist rate. Although CS approaches have been proposed for ultrasound (US) imaging with promising results, practical implementations are hard to achieve due to the impossibility to mimic random sampling on a US probe and to the high memory requirements of the measurement model. In this paper, we propose a CS framework for US imaging based on an easily implementable acquisition scheme and on a delay-and-sum measurement model.

We propose a Gauss-Legendre quadrature based sampling on the rotation group for the representation of a band-limited signal such that the Fourier transform (FT) of a signal can be exactly computed from its samples. Our figure of merit is the sampling efficiency, which is defined as a ratio of the degrees of freedom required to represent a band-limited signal in harmonic domain to the number of samples required to accurately compute the FT. The proposed sampling scheme is asymptotically as efficient as the most efficient scheme developed very recently. For the computation of FT and inverse FT, we also develop fast algorithms of complexity similar to the complexity attained by the fast algorithms for the existing sampling schemes. The developed algorithms are stable, accurate and do not have any pre-computation requirements. We also analyse the computation time and numerical accuracy of the proposed algorithms and show, through numerical experiments, that the proposed Fourier transforms are accurate with errors on the order of numerical precision.

This paper reviews and extends our previous work to enable fast axonal diameter mapping from diffusion MRI data in the presence of multiple fibre populations within a voxel. Most of the existing mi-crostructure imaging techniques use non-linear algorithms to fit their data models and consequently, they are computationally expensive and usually slow. Moreover, most of them assume a single axon orientation while numerous regions of the brain actually present more complex configurations, e.g. fiber crossing. We present a flexible framework, based on convex optimisation, that enables fast and accurate reconstructions of the microstructure organisation, not limited to areas where the white matter is coherently oriented. We show through numerical simulations the ability of our method to correctly estimate the microstructure features (mean axon diameter and intra-cellular volume fraction) in crossing regions.

We develop a novel sampling theorem for functions defined on the three-dimensional rotation group SO(3) by connecting the rotation group to the three-torus through a periodic extension. Our sampling theorem requires 4L3samples to capture all of the information content of a signal band-limited at L, reducing the number of required samples by a factor of two compared to other equiangular sampling theorems. We present fast algorithms to compute the associated Fourier transform on the rotation group, the so-called Wigner transform, which scale as O(L4), compared to the naive scaling of O(L6). For the common case of a low directional band-limit N, complexity is reduced to O(NL3). Our fast algorithms will be of direct use in speeding up the computation of directional wavelet transforms on the sphere. We make our S03 code implementing these algorithms publicly available.

In recent works, compressed sensing and convex optimization techniques have been applied to radio-interferometric imaging showing the potential to outperform state-of-the-art imaging algorithms in the field. We review our latest contributions, which leverage the versatility of convex optimization to both handle realistic continuous visibilities and offer a highly parallelizable structure paving the way to significant acceleration of the reconstruction and high-dimensional data scalability. The new algorithmic structure, promoted in a new software PURIFY (beta version), relies on the simultaneous-direction method of multipliers (SDMM). The performance of various sparsity priors is evaluated through simulations in the continuous visibility setting, confirming the superiority of our recent average sparsity approach SARA.

Inspired by the recently proposed Magnetic Resonance Fingerprinting technique, we develop a principled compressed sensing framework for quantitative MRI. The three key components are: a random pulse excitation sequence following the MRF technique; a random EPI subsampling strategy and an iterative projection algorithm that imposes consistency with the Bloch equations. We show that, as long as the excitation sequence possesses an appropriate form of persistent excitation, we are able to achieve accurate recovery of the proton density, T1, T2and off-resonance maps simultaneously from a limited number of samples.

Image recovery in optical interferometry is an ill-posed nonlinear inverse problem arising from incomplete power spectrum and bispectrum measurements. We review our previous work, which reformulates this nonlinear problem in the framework of tensor recovery and studies two different approaches to solve it: one is nonlinear and nonconvex while the other is linear and convex. We extend the linear convex procedure to account for signal sparsity and we also present numerical simulations that show the improvement in the quality of reconstruction of sparse images when including a sparsity prior.

Inspired by the recently proposed magnetic resonance fingerprinting (MRF) technique, we develop a principled compressed sensing framework for quantitative MRI. The three key components are a random pulse excitation sequence following the MRF technique, a random EPI subsampling strategy, and an iterative projection algorithm that imposes consistency with the Bloch equations. We show that, theoretically, as long as the excitation sequence possesses an appropriate form of persistent excitation, we are able to accurately recover the proton density, T1, T2, and off-resonance maps simultaneously from a limited number of samples. These results are further supported through extensive simulations using a brain phantom.

Validation is arguably the bottleneck in the diffusion magnetic resonance imaging (MRI) community. This paper evaluates and compares 20 algorithms for recovering the local intra-voxel fiber structure from diffusion MRI data and is based on the results of the "HARDI reconstruction challenge" organized in the context of the "ISBI 2012" conference. Evaluated methods encompass a mixture of classical techniques well known in the literature such as diffusion tensor, Q-Ball and diffusion spectrum imaging, algorithms inspired by the recent theory of compressed sensing and also brand new approaches proposed for the first time at this contest. To quantitatively compare the methods under controlled conditions, two datasets with known ground-truth were synthetically generated and two main criteria were used to evaluate the quality of the reconstructions in every voxel: correct assessment of the number of fiber populations and angular accuracy in their orientation. This comparative study investigates the behavior of every algorithm with varying experimental conditions and highlights strengths and weaknesses of each approach. This information can be useful not only for enhancing current algorithms and develop the next generation of reconstruction methods, but also to assist physicians in the choice of the most adequate technique for their studies.

We discuss a novel sparsity prior for compressive imaging in the context of the theory of compressed sensing with coherent redundant dictionaries, based on the observation that natural images exhibit strong average sparsity over multiple coherent frames. We test our prior and the associated algorithm, based on an analysis reweighted formulation, through extensive numerical simulations on natural images for spread spectrum and random Gaussian acquisition schemes. Our results show that average sparsity outperforms state-of-the-art priors that promote sparsity in a single orthonormal basis or redundant frame, or that promote gradient sparsity. Code and test data are available at https://github.com/basp-group/sopt.

We propose a novel algorithm for source signals estimation from an underdetermined convolutive mixture assuming known mixing filters. Most of the state-of-the-art methods are dealing with anechoic or short reverberant mixture, assuming a synthesis sparse prior in the time-frequency domain and a narrowband approximation of the convolutive mixing process. In this paper, we address the source estimation of convolutive mixtures with a new algorithm based on i) an analysis sparse prior, ii) a reweighting scheme so as to increase the sparsity, iii) a wideband data-fidelity term in a constrained form. We show, through theoretical discussions and simulations, that this algorithm is particularly well suited for source separation of realistic reverberation mixtures. Particularly, the proposed algorithm outperforms state-of-the-art methods on reverberant mixtures of audio sources by more than 2 dB of signal-to-distortion ratio on the BSS Oracle dataset.

We study the impact of sampling theorems on the fidelity of sparse image reconstruction on the sphere. We discuss how a reduction in the number of samples required to represent all information content of a band-limited signal acts to improve the fidelity of sparse image reconstruction, through both the dimensionality and sparsity of signals. To demonstrate this result, we consider a simple inpainting problem on the sphere and consider images sparse in the magnitude of their gradient. We develop a framework for total variation inpainting on the sphere, including fast methods to render the inpainting problem computationally feasible at high resolution. Recently a new sampling theorem on the sphere was developed, reducing the required number of samples by a factor of two for equiangular sampling schemes. Through numerical simulations, we verify the enhanced fidelity of sparse image reconstruction due to the more efficient sampling of the sphere provided by the new sampling theorem.

We propose a novel compressed sensing technique to accelerate the magnetic resonance imaging (MRI) acquisition process. The method, coined spread spectrum MRI or simply s(2)MRI, consists of premodulating the signal of interest by a linear chirp before random k-space under-sampling, and then reconstructing the signal with nonlinear algorithms that promote sparsity. The effectiveness of the procedure is theoretically underpinned by the optimization of the coherence between the sparsity and sensing bases. The proposed technique is thoroughly studied by means of numerical simulations, as well as phantom and in vivo experiments on a 7T scanner. Our results suggest that s(2)MRI performs better than state-of-the-art variable density k-space under-sampling approaches.

Radio interferometry is a powerful technique for astronomical imaging. The theory of compressed sensing (CS) has been applied recently to the ill-posed inverse problem of recovering images from the measurements taken by radio interferometric telescopes. We review novel CS radio interferometric imaging techniques, both at the level of acquisition and reconstruction, and discuss their superior performance relative to traditional approaches. In order to remain as close to the theory of CS as possible, these techniques necessarily consider idealised interferometric configurations. To realise the enhancement in quality provided by these novel techniques on real radio interferometric observations, their extension to realistic interferometric configurations is now of considerable importance. We also chart the future direction of research required to achieve this goal.

Incoherence between sparsity basis and sensing basis is an essential concept for compressive sampling. In this context, we advocate a coherence-driven optimization procedure for variable density sampling. The associated minimization problem is solved by use of convex optimization algorithms. We also propose a refinement of our technique when prior information is available on the signal support in the sparsity basis. The effectiveness of the method is confirmed by numerical experiments. Our results also provide a theoretical underpinning to state-of-the-art variable density Fourier sampling procedures used in MRI.

We consider images probed through incomplete and noisy Fourier coverages, both in the context of radio interferometry (RI) and of magnetic resonance imaging (MRI). We show that the quality of signal reconstruction can be significantly enhanced by the introduction of a linear chirp modulation, which induces a spread spectrum phenomenon.