Convergent Regularization in Inverse Problems and Linear Plug-and-Play Denoisers

Foundations of Computational Mathematics https://doi.org/10.1007/s10208-024-09654-x Convergent Regularization in Inverse Problems and Linear Plug-and-Play Denoisers Andreas Hauptmann1,2 · Subhadip Mukherjee3 · Carola-Bibiane Schönlieb4 · Ferdia Sherry4 Received: 14 December 2023 / Revised: 25 March 2024 / Accepted: 25 March 2024 © The Author(s) 2024 Abstract Regularization is necessary when solving inverse problems to ensure the wellposedness of the solution map. Additionally, it is desired that the chosen regularization strategy is convergent in the sense that the solution map converges to a solution of the noise-free operator equation. This provides an important guarantee that stable solutions can be computed for all noise levels and that solutions satisfy the operator equation in the limit of vanishing noise. In recent years, reconstructions in inverse problems are increasingly approached from a data-driven perspective. Despite empirical success, the majority of data-driven approaches do not provide a convergent regularization strategy. One such popular example is given by iterative plug-and-play (PnP) denoising using off-the-shelf image denoisers. These usually provide only convergence of the PnP iterates to a fixed point, under suitable regularity assumptions on the denoiser, rather than convergence of the method as a regularization technique, thatis under vanCommunicated by Teresa Krick and Hans Munthe-Kaas. Invited paper associated to the FoCM 2021 Online Seminar lecture Machine Learned Regularization for Solving Inverse Problems presented by Carola-Bibiane Schönlieb in April 2021. B Carola-Bibiane Schönlieb Andreas Hauptmann Subhadip Mukherjee Ferdia Sherry 1 Research Unit of Mathematical Sciences, University of Oulu, Oulu, Finland 2 Department of Computer Science, University College London, London, UK 3 Department of Electronics and Electrical Communication Engineering, Indian Institute of Technology Kharagpur, Kharagpur, India 4 Department of Applied Mathematics and Theoretical Physics, University of Cambridge, Cambridge, UK 123 Foundations of Computational Mathematics ishing noise and regularization strength. This paper serves two purposes: first, we provide an overview of the classical regularization theory in inverse problems and survey a few notable recent data-driven methods that are provably convergent regularization schemes. We then continue to discuss PnP algorithms and their established convergence guarantees. Subsequently, we consider PnP algorithms with learned linear denoisers and propose a novel spectral filtering technique of the denoiser to control the strength of regularization. Further, by relating the implicit regularization of the denoiser to an explicit regularization functional, we are the first to rigorously show that PnP with a learned linear denoiser leads to a convergent regularization scheme. The theoretical analysis is corroborated by numerical experiments for the classical inverse problem of tomographic image reconstruction. Keywords Inverse problems · Variational regularization · Data-driven learning · Plug-and-play denoising Mathematics Subject Classification 47A52 · 46N10 · 65F22 1 Introduction Inverse problems deal with the estimation of an unknown model parameter x ∗ ∈ X from its noisy and indirect measurement y δ ∈ Y given by y δ = Ax ∗ + e. (1) We consider the case where X and Y are (potentially infinite dimensional) separable Hilbert spaces and A : X → Y is a bounded linear operator. X and Y are endowed with inner products ·, · X and ·, ·Y , inducing the norms · X and ·Y , respectively. The measurement noise level is bounded by δ, i.e., eY ≤ δ. The clean measurement is denoted by y 0 . The inverse problem in (1) is considered ill-posed in the sense of Hadamard, if either injectivity or surjectivity of the forward operator, or stability of the solution map is violated. For instance, if A is a compact operator with an infinite-dimensional range, then surjectivity and stability are not satisfied. This is, for example, the case for the ray transform operator that underlies many applications in medical imaging, such as computed tomography (CT) and positron emission tomography (PET) [35, 36]. The study of inverse problems usually assumes ill-posedness, as we will also do in the following. To address ill-posedness, one needs to introduce a general concept for stable and unique solvability for an inverse problem of the form (1). Due to the aforementioned ill-posedness, we can not guarantee the recovery of the true solution x ∗ for all measurements and hence we first need the concept of a generalized solution. A common approach is to search for solutions that are closest to the measured data with respect to a suitable data discrepancy term f : Y × Y → R+ , such as the (squared) distance 123 Foundations of Computational Mathematics in the norm, i.e., f (Ax, y δ ) = Ax − y δ 2Y . Then we search for  x ∈ X such that f (A x , y δ ) ≤ f (Ax, y δ ) for all x ∈ X . (2) (2) implies that  x is closest to the measured data with respect to f , which deals with the violation of surjectivity by disregarding components of y δ in the co-kernel of A. Furthermore, if A has a non-trivial null space, then  x is not unique. To obtain a unique solution, one can define the minimum norm solution as x † = arg min{x X : x minimizes f (Ax, y δ )}. (3) x∈X The element x † can now be considered a desirable generalized solution to (1). When f and  ·  X are given by the squared L 2 -norm, we call x † the least-squares minimumnorm solution and can define a mapping A† : Y → X , such that x † = A† y δ . In fact, the mapping A† defines what is referred to as the Moore-Penrose pseudo-inverse. Unfortunately, if the operator A is compact, then A† will be unbounded and as such does not take care of the stability problem in the presence of noise in the data. This is where the concept of regularization becomes important, as we will discuss next. Regularization theory considers specifically designed solution maps to deal with the stability issue. Such a solution map R(·; λ) : Y → X , also called a reconstruction operator, is expressed as a parametric map that produces a solution estimate of x ∗ given y δ . Here, the parameter λ depends on the noise level δ and the measured data y δ , which we denote explicitly by the mapping λ = λ(δ, y δ ). In this paper, we are specifically interested in the notion of convergent regularization which can be understood as convergence of the reconstruction operator when the noise level δ tends to zero. More specifically, we want that when the noise level δ → 0, then λ(δ, y δ ) → λ0 ≥ 0, and the reconstruction operator R(y δ ; λ) converges to a generalized solution of the noiseless operator equation Ax = y 0 . (4) A family of such reconstruction operators R(y δ ; λ) can be formulated in the framework of variational regularization (see Sect. 2.1.3 for more details) by defining them as th (...truncated)