NPN: Non-Linear Projections of the Null-Space for Imaging Inverse Problems
Authors: Roman Jacome
Romario Gualdrón-Hurtado
León Suárez-Rodríguez
Henry Arguello
Venue: NeurIPS 2025 poster ()
Keywords: Inverse problems
machine learning
computational imaging
optimization.
Primary area: optimization
TL;DR: Our work introduces a solver-agnostic null-space regularizer with gains across inverse problems.

Abstract: Imaging inverse problems aim to recover high-dimensional signals from undersampled, noisy measurements, a fundamentally ill-posed task with infinite solutions in the null-space of the sensing operator. To resolve this ambiguity, prior information is typically incorporated through handcrafted regularizers or learned models that constrain the solution space. However, these priors typically ignore the task-specific structure of that null-space.  In this work, we propose Non-Linear Projections of the Null-Space (NPN), a novel class of regularization that, instead of enforcing structural constraints in the image domain, promotes solutions that lie in a low-dimensional projection of the sensing matrix's null-space with a neural network.
Our approach has two key advantages: (1) Interpretability: by focusing on the structure of the null-space, we design sensing-matrix-specific priors that capture information orthogonal to the signal components that are fundamentally blind to the sensing process. (2) Flexibility: NPN is adaptable to various inverse problems, compatible with existing reconstruction frameworks, and complementary to conventional image-domain priors. We provide theoretical guarantees on convergence and reconstruction accuracy when used within plug-and-play methods. Empirical results across diverse sensing matrices demonstrate that NPN priors consistently enhance reconstruction fidelity in various imaging inverse problems, such as compressive sensing, deblurring, super-resolution, computed tomography, and magnetic resonance imaging, with plug-and-play methods,  unrolling networks, deep image prior, and diffusion models.

NPN: Non-Linear Projections of the Null-Space for Imaging Inverse Problems
Roman Jacome Romario Gualdrón-Furtadoo; Leon Suarez', Henry Arguello
DDeparmeen of Electrical, Electronics, and Telecommunications Engineering
*Department of Systems Engineering and Informatics
Universidad Industrial de Santander, Colombia. 680002
ead...dd..com co, henarfuduis edu..co
Abstract
Imaging inverse problems aim to recover high-dimensional signals from undersam- pled, noisy measurements, fundamentally ill- ill-posed task with infinite solutions in the null-space of the sensing operator. To resolve this ambiguity, prior in- formation is typically incorporated through handerafted regularizers learned models that constrain the solution space. However, these priors typically ignore the task-specific structure of that null-space. In this work, propose Non-Linear Projections of (he Null-Space (NPN), novel class of regularization that, instead of enforcing structural constraints in the image domain, promotes solutions that lie in a oow-dimesioonal projection of the sensing matrix's 's null-space with a neural network. Our approach has key advantages: (1) Interpretability: by focusing on the structure of the null-space, we design sensing- sensing-matrix mntrix--peecific priors that capture information orthogonal to the signal components that are fundamentally blind to the sensing process. (2) Flexibility: NPN is adaptable to various inverse problems, compatible with existing reconstruction frameworks, and complemen- tary to conventional image-domain priors. We provide theoretical guarantees on convergence and reconstruction accuracy when used within plug-and-play methods. Empirical results across diverse sensing matrices demonstrate that NPN priors consistently enhance reconstruction fidelity various imaging inverse problems, such as compressive sensing, deblurring, super-resolution. computed tomography, and magnetic resonance imaging, with plug-and-pla methods, unrolling networks, deep image prior, and diffusion models.
1 Introduction
Inverse problems involve reconstructing an unknown signal from noisy, corrupted, or undersampled observations, making the recovery process generally non-invertible and ill-posed. This work focuses on linear inverse problems of the form = Hxx + w E Rm, where x* E R is the target high- dimensional signal, H E Rmxn is the sensing matrix (with m & n2), y I" represents the low-dimensional measurements, and w NO, I² is additive Gaussian noise. Numerous imaging tasks rely on this principle, including image restoration-sul as deblurring, denoising, inpainting, and super-resolution (SR) [21] (structured Toeplitz sensing matrices) well as compressed sensing (CS) [61, 4] (dense sensing matrices) and medical imaging applications like magnetic resonance imaging (MRI) [36] (undersampled Fourier matrices) matrices). The challenge on the recovery task lies in the ill ill-posed nature of the inverse problem due to the non-trivial null-space of the sensing matrix H leading infinite solutions. Therefore, there a need to incorporate a signal prior to the reconstruction framework. Under this idea, variational approaches formulate the signal estimator
computed tomography (CT) [59] (Radon
X arg min g(x) + Xh(x)
(1)
x
where g(-) is the data fidelity term and f(-)) is regularization function based on some prior of x".
'Equal contribution.
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

(21g)
y only Popttnd
10
Tea
1.0g0 60
-15
5.0
60 15 0.0
Figure I: Geometric comparison of subspace-prion learning versus direct reconstruction in TR toy example. (a) In the low- dimensionn projection space, the learned mapping G' (y) trained on points inside the unit circle, closely matches the true null-space projection Sx' for both training (solid) and test (eemi-transparent inputs, whereas the reconstriccnnt estimate Ka projected into is significantly inaccurate. (b) In the original signal domain, the measurements Hx" = define two intersecting planes.
One of the most common priors in imaging inverse problems is sparsity, which assumes that images are compressible in a given basis [4]. Other widely used priors include low-rank structures [17] and smoothness priors [19]. Additionally, plug-and-play (PnP) priors [55], which traces back its roots from proximal algorithms [40], where these operators, usually defined by analytical models of the underlying signals [61], are replaced by a general-purpose purpose image denoiser [5 J. This approach allows the integration of classical image denoiser [37] such as BM3D [14], NLM [2), RF [16] and current deep learning (DL) denoisers [62]. The idea behind DL-based denoisers is to train deep neural network (DNN) that maps from the noisy observation to the clean image [3. 53, 52, 24, 23]. Another (earning-based approaches are based on the null-space of H by embedding the sensing operator's 's structure directly into learned networks. In particular, Null-Space Networks exploit the decomposition of signal into measurement and null-space components, learning a corrective mapping over all null- space modes to enhance interpretability and accuracy [45]. To improve robustness to measurement noise, [6] introduced separate range-space and null-space networks that denoise both components before recombination. Variants of this range-null decomposition have been applied in diffusion- based restoration [11, 57, 58], GAN-prior methods [56], uggorithm-unrolling architectures [5], and self-supervised schemes [7, 10], consistently leveraging the full null-space projector to achieve high-fidelity reconstructions.
However, existing learned priors typically promote reconstructions that lie within the subspace spanned by clean training data, without explicitly accounting for the null-space of the sensing matrix H. While the data fidelity term g(x) enforces consistency with the measurements, it does not sufficiently constrain the null-space components of the solution, especially in the presence of noise, often resulting in suboptimal reconstructions. In this work, we introduce a novel class of regularization, termed Non-Linear Projections of the Null-Spuce (NPN), which directly promotes solutions within a low-dimensional subspace of the null-space of H that is, within the space of vectors orthogonal to the rows of H. Our method identifies compact subspace of the null-space by selecting only the most informative directions and trains a neural network to predict their coefficients directly from measurements. By restricting corrections to this learned subspace, we concentrate regularization on unobserved features most predictive of the true signal. This subspace plays a critical role in addressing the ill- posedeesssss of inverse problems, where conventional methods often struggle due to the lack of constraints in directions invisible to the measurement operator. To enable this, we design a projection matrix S whose rows lie in the null-space of H. constructed using either orthogonalization techniques or analytical designs depending on the structure of H. A neural network is trained to estimate this null-space projection from the measurements y. providing a non-linear prior that is both data-adaptive and model-aware. We further propose a joint optimization framework in which both the projection matrix S and the network are learned simultaneously, allowing the projection matrix to adapt during training in task-specific manner.
Our approach offers two key advantages for imaging inverse problems: interpretability and fiexibility. Interpretability: By leveraging a non-linear neural network G* to learn projection onto a low- dimensional subspace within the null-space of H. we focus on recovering signal components that lie outside the range of the sensing matrix. Learning this projection, rather than directly estimating the full signal x", reduces complexity through dimensionality reduction while maintaining clear connection to the geometry of the inverse problem. Flexibility: the learned prior can be readily incorporated into a wide variety of reconstruction algorithms and image priors that can be adapted to other imaging inverse problems. To illustrate the interpretability advantage, Figure ((a) shows
2

---

that the learned non-linear estimator G closely matches the true null-space projection Sx* (green surface), while the projection of the direct reconstruction X, onto the same subspace spanned by S results in significant errors. Moreover, when evaluating grid in the range [2, 4J2, our method experiences only a minor increase in estimation error, while the direct reconstruction approach drastically amplifies these errors. In Figure ((b), we demonstrate that integrating our non-linear subspace prior into the inverse reconstruction process effectively regularizes the solution space. The learned subspace helps uniquely resolve the inverse problem, providing a solution close to the true signal x*, A zoomed zone near the true solution is shown, where there is small estimation eтrоr with respect to the true value; this is due the inherent network error in the true subspace. Note that here we did not use any prior hxx), which can reduce the estimation error. See the Appendix A.1 for more details on the setting obtain this example. We develop theoretical analysis showing benefits in the convergence rate when the NPN regularizer is used in PnP algorithms. The theory states that the algorithms have significant speed-up (with respect to the non-regularizer counterpart) in convergence improvement zone (CIZ), which takes into account the inherent estimation error of the learned prior. Additionally, we showed that the NPN regularizer is minimized up to a constant depending on the error of the learned prior if the algorithm reaches optimum values. The theory is validated in a comprehensive evaluation of the method in five imaging inverse problems: CS, MRI, deblurring, SR, and CT. Although our theoretical findings are based on PnP methods, we also validated the NPN regularization in unrolling networks, deep image prior, and diffusion models.
samples within a 2 x 2
2 Related work
2.1 Variational reconstructions methods
Variational methods solve (1) typically via proximal-gradies schemes that alternate gradient step on the data fidelity term, g(x) = |H y1l.. with a proximal step for h. A classical choice is Tikhonov regularization for smoothness, h(x) |[x 11.. where L is derivative operator identity [19]. Alternatively, sparsity priors use an é li-penalty in transform domain (e.g. wavelets) to promote compressibility of x [4]. Algorithms such as ISTA and its accelerated variant FISTA have been widely used to solve é--penalized reconstructions with provably faster convergence rates [1]. More recently, PnP replaces the proximal operator of h with a generic denoiser x* = D, (x aH (Hx* y)). thus, leveraging sophisticated image priors without an explicit analytic penalty [55, 62]. PnP with both classical denoisers and deep-learning models achieves state-of-the-an results, yet it leaves the null-space of H uncontrolled: any component in the null-space of H may be arbitrarily modified by the denoiser.
to penalize energy and ensure well-posed
2.2 Null-space learned reconstruction
Harnessing the sensing model into a learning- reconstruction network enables more accurate reconstructions [45]. Particularly, null-space networks (NSN) [45] harness the range-null-space decomposition (RNSD), which states that a vector x E R" is decomposed as x = (x(x) + (x) where (x) = H' with H denoting the Moore-
pseudoinverse, is the projection onto the range space of H and P, (x) = (I H H)) is the null-space projection operator. NSN exploits this structure by learning neural network R R + IR" such that the reconstruction becomes H + (I- H'')RRH''S (H) This approach benefits the interpretability of the reconstruction. However, this method does not take into account the inherent noise of the measurements; thus, [6] introduces deep decomposition networks (DDN), a more robust formulation is presented where two models R, and R, are employed to improve recovery performance. Consider two architectures, DDN- independent (DDN-I) & H'y + (, (H))) + P,, (R,, (H'Y)) and DDN-Cascade (DDN-C) H'y + Pr (R, (H' y)) + P,, (R,, (H, P. (R, (H y)))). The range-null-space decomposition has also been used to enhance data fidelity in diffusion models for image restoration [1 111, 57, 58), generative adversarial network priors [56], unfolding networks [5], and self supervised learning [7, 10]. In these works, the sensing matrix H structure is incorporated into the reconstruction network through a learned-based RNSD that allows high-fidelity reconstructions. Different from these approaches, our method does not apply the full null-space projection operator. Instead, we first identify a compact subspace of Null(H) by selecting a projection matrix S E RR whose rows span the most informative null-space directions. We then train network G* 1 + R" to predict the
3

---

coefficients y, = Sx directly from the measurements. On the other hand, the learned NSN [45] is optimized to improve the recovery performance of determined regularized inverse problem (i.e.. Tikhonov-based solutions), which makes it algorithm-specife and does not work as plug-in for other recovery methods. In our case, since we optimize the network offline, only with the knowledge of H and S, it can be easily integrated in wide range of recovery methods.
3 Method
In our approach, we design learned priors promoting solutions in a low-dimensional subspace of the null-space of H. First, let's define some useful properties.
Definition (Null-Space). The null-space of a matrix H is defined as
Hu)l(H) = (x E R" Hx = 0} = fx x 1 hy. Vj E (1....,m,}
Thus, we consider projection matrix S R2xx" with p S (n-m), with rows orthogonal H rows, implying that s; 1 h, Vi E (1.... p). .). Vjj E {1..... m.. Based on Def. 1, s; E Null(H) meaning that any projection y, Sx lies onto a low-dimensional subspace of the null-space of H. Based on this observation, we propose to learn a data-driven prior G(-) restricted to the low-dimensional null-space of H, such that G(y) = Sx*. Specifically, we select a projection matrix S E Rexa whose rows span a subspace of NH). Consequently, we solve
x arg min g(x) + Ah(x) 7 ||G (y) SX
(2)
and G* 1Rm + IRP is neural network trained map the measurements y into the low-dimensional subspace Sx*, 7 is a regularization parameter, the regularizer 6 aims to promote solutions on the learned manifold in the null-space of H. Our framework introduces a novel regularization strategy that embeds data-driven models into inverse- inverse-problem solvers by constraining solutions to the nonlinear oiw-dimensional manifold induced by G in contrast to existing learned priors that restrict reconstructions to the range of a pre- pre-trained restoration or denoising network [52, 24, 291. Note that we used the Fuclidean norm in the regularizer; however, since the network G has some error with respect to Sx", one could use more robust function such as the Huber loss or weighted norm. Nevertheless, in our experiments, the Euclidean norm works well by adjusting the hyperparameter y. One interpretation of the proposed regularization is that it improves the data-fidelic term ((x) as it promotes fow-dimensional projections of blind signal features to H. Our approach is closely related to NSN-based methods [6, 45, 58], as those methods aim to regularize deep learning-base recovery networks, harnessing the RNSD. our case, making analogy with these models, we can view the reconstruction in (2) as x= H'y + Sy) prior,,, where the projection onto the null-space (-) is replaced by the range-space S which promotes solution lying in the S (n m) most informative null-space modes instead of the entire null-space operator.
3.1 Design of the matrix
To design the matrix S, it is necessary to analyze the structure of H depending on the inverse problem. Below are insights based on the sensing matrix structure for exploiting the null-space in our prior.
Compressed Sensing (CS): In CS, the sensing matrix H E Rm is typically dense and randomly generated. Previous approaches [20, 38] use the remaining (n m) rows of full-rank Hadamard, Gaussian, or Bernoulli matrix as S. Due to the lack of inherent structure in such matrices, we adopt an orthogonalization strategy based on the classical QR decomposition for designing the matrix (see Alg. 1 in Appendix A.2).
Magnetic Resonance Imaging (MRI): In MRI, the forward operator H corresponds to a discrete 2D discrete Fourier transform (DFT) undersampled, where only a subset of frequency components (k-space lines) is acquired during the scan. Specifically, let J {f, ..... f denote the full set of orthonormal rows of the 2D DFT matrix. The sensing matrix H then consists of m < n selected rows, i.e., (h, ...., h == F, F. These rows define the measurements taken in the Fourier domain. To construct the null-space projection matrix S, we exploit the fact that the remaining rows in F T-those not used in H span the null-space of the sampling operator. Thus, we define {s] ,...., 8, = Fÿ. where Fí = F \ Fú is the complement of the sampled frequencies. Because the DFT matrix is orthonormal, these complementary rows are guaranteed to be orthogonal to the measurement space and form natural basis for the null-space of H.

---

Computed Tomography (CT); In parallel beam limited- angle CT, the forward operator samples the Radon transform only at a subset of projection angles. Let be the full discrete angle set, O₂ C ④ the acquired angles, R is the discrete Radon transform matrix with rows indexed by θ and H Pe, R the forward operator; we define S directly as the complement of the acquired angles, S = Peg R. with 6 = θ b so S stacks the rows of R corresponding to the non-acquired angles.
Structured Toeplitz matrices (Deblurring and Super-Resolution. The forward model H is built upon a Toeplitz matrix based on the convolution kemell denoted as H(ii, i j] = híj] with i= 1. m and j I......2 n. The action of H corresponds a linear filtering process, attenuating high-frequency components. From frequency point of view, the matrix S should block the low frequencies sampled by H. Thus, we design the matrix as Sii, i j] and j L..... .R.. In super-resolution. the sensing matrix is H = DB where B E R"** is a convolution matrix build with a low-pass filter b and D TR m is decimation matrix denoting V V the super-resolution factor (SRF). We construct the matrix S similarly to the deblurring case, i.e., S 2 + jl = 1 - 1-b[j] with i = I...., m and =1, 1....,
h[j] with i 1..... m
3.2 Learning the NPN Prior
Given the design of the projection matrix S, which spans a structured low-dimensional subspace orthogonal to the measurement operator H. we train the network G to estimate the null-space component x" of the signal. To further improve the representation power of the NPN prior, we propose to jointly optimize the neural network G and the projection matrix S. While the initial design of ensures that its rows lie in the null-space of the measurement operator H. optimizing S during training allows the model to discover task-adaptive subspace that best complements the measurements. This formulation still preserves the orthogonality between H and S, while improving the quality of the previously learned. The optimization objective is
G. arg min Ex, y ||G(Hxx) - Sxx +1llx - A'Ax A - -1 G
(3)
where T X1, A2 > 0 control the trade-off between estimation accuracy and orthog- onality enforcement. The first term trains G to predict the projection Sx from the measurements y, following the MMSE objective. The second term enforces near-orthogonality between the row spaces of S and H by pushing S S + H H se I. The third term promotes full-rank behavior and numerical stability of the combined system matrix A, preventing collapse or redundancy in the learned subspace. Importantly, this formulation enables end-to-end learning of a null-space mull-space-awan regularizer, where the matrix S is initialized using principled designs (e.g.. QR orthogonalization or frequency complements in MRI), but is then refined during training to maximize consistency with the true signal statistics and the learned estimator G. This formulation has great benefits in cases of non-structured dense matrices, such as those in CS; however, for well-defined matrices, such as Fourier-based or Toeplitz matrices, finding orthogonal complements has a straightforward analytical solution, and it is not required to jointly optimize it.
Theoretical analysis
We analyze how the convergence property of a PnP algorithm is affected under this regularization. Without loss of generalization, we focus on one of the most common PnP approaches, which is based on proximal gradient methods. The iterations are given by
T(x() = Dz (x a (Vg(X") Vo(x))))
where (.) is the denoiser operator and a is an hyperparameter modeling the noise variance and L = 1, ... L with maximum number of iterations L.
Based on this formulation, we analyzed two main aspects: and (ii) the convergence behavior of the regularization function Φ(x*+1) = [[G'(y) -x - Our theoretical developments leverage the restricted isometry property defined over specific Riemannian manifold MD, induced by the image-space of the denoiser D,. For the denoiser to exhibit isometric- like properties, certain criteria must be satisfied, including boundedness, Lipschitz continuity, and a low-rank Jacobian. These properties can be guaranteed through spectral normalization during denoiser training [43]. Such assumptions are commonly employed in the convergence analysis of
(i) the convergence rate of the algorithm.
5

---

iterative projection methods [46] and have been suitably adapted for PnP convergence analyses [31]. Additionally, we introduce assumptions regarding the estimation error of the model G*, assuming a Gaussian error distribution. Furthermore, we define CI7. in the algorithm iterations based on the estimation error norm, For the guarantees, we consider the noiseless case, i.e., w O.
Definition 2 (Restricted Isometry Property [46]). Lei MD C. RR" be Riemannian manifold given by the denoiser's 's image space Imn(Dr), ). thus TRpx0 satisfies the restricted isometry property with respect to MD with a restricted isometry constant (RIC) 2M E (0,, 1) if for all x, n. E MD with
(1 - ) 21: S |S(x 2 S (1 + A ) z 2
(5)
Assumption 1 (Prior mismatch). The trained model G using (3), we assume that
frath vulur
G G' Sx + +N(Hx")
(6)
Nov-lineer crent term
Thus, considering that the nonlinear operator N is K K-Lipschitz continuous, and H satisfies Definition 2 with a constant A thus we have (N( )11 S K(1 + A 1|
Definition 3 (Convergence Improvement Zone (CIZ)) by d(x²)), We define a zone where the proposed NPN prior improves the convergence of the PnP algorithm. For that, S satisfies the RIP with RIC AMD 1), the CIZ by & are the iterations Lo {I,...; Lo} where Log S L such that the network estimation error ((Hx) for all E Lus satisfies that
[[N(Hx) 11 < |SX* Sx'1 = (S - x') 112 S (1+ (L+AM)))IE* -x'11
Assumption 2 (Bounded denoiser). We consider that for x, 2 E R", Do is a bounded denoiser, with
constant 8 > if
(x) Da(z)]1 & (1 6][[ 2
We are now equipped to develop the first theoretical benefit of the proposed method.
Theorem (PnP-NPN Convergence). Consider the fidelity term g(x) - lly H and assume the denoiser Dr satisfies Assumption 2. Let the matrix S be constructed according to (3) and satisfy the RIP condition (Definition 2) with constant MD E [0, 1). Then, for a finite number of iterations
{= (=I,...; Ló within the CIZ, the residual ||x* ++]
xx1 x212 decays linearly with rate
p _ (1 + 8) (III [ S S))l2 + (1+AA <1
The proof of the theorem can be found in Appendix A.3. A key insight from Theorem 1 is the role of the CIZ. L.. which characterizes the subset of iterations where the proposed regularizer outperforms provides improved convergence. Specifically, due to the inherent mismatch between the ground-truth projection Sx" and its learned estimate GO (Hx"), ), the NPN prior only effective while the projected estimate Sx remains closer to Sx' than the residual error N N(Hx") ). Outside this zone, the regularizer may no longer provide beneficial guidance to the reconstruction. Nevertheless, thanks to the design of the matrix S either through analytical construction or data-driven optimization, is guaranteed to be orthogonal to the rows of the sensing matrix H. This orthogonality ensures that the operator
norm ||I - (H H + S) remains small for an appropriate step size a. Furthermore, when the spectral norm of S and its restricted isometry constant A S Mp are sufficiently low, the acceleration factor p falls below one. This guarantees that the PnP-NPN algorithm will converge a fixed point within the zone Lg. thereby validating the theoretical and practical benefits of incorporating the NPN prior into the reconstruction framework.
The second analysis of the proposed approach is the convergence of the regularization. Theorem 2 (Convergence of NPN Regularization), Consider the tterations of the PnP-NPN algorithm defined in for i 1..... L. Assume that the estimation error of the trained network G satisfies Assumption I, and that both S and H. satisfy the Restricted Isometry Property (RIP) over the manifold
MD with constants A A AHP
MD E [0, 1). Further assume that the residual term N(Hx*) is K-
Lipschitz continuous, thus [[N(Hx") )13 < K(1 + A Mp >||x* 112. Then, after E iterations, the NPN
regularization term satisfies the following upper bound:
[[G(y) - -x 112 < C1 ||x" - x x113 + X
21
- + C21lx -x 1
(8)
where C1 = ( 1 + K A ||x (13) (1+ AM.)
C2 (1+ (LD) )² (1+

---

The proof of Theorem 2 is provided in the technical Appendix A.4. This result shows that the regularization value ((x*+1) decreases as the reconstruction error ||x* - -x112 diminishes over the course of the iterations. From Theorem we know that ||x* x [2 is reduced, thereby ensuring a monotonic decrease in the regularization term. Moreover, since Theorem 1 guarantees convergence of the sequence {x'} a fixed point, the difference ||x* x asymptotically vanishes as e + 0C. Therefore, in the limit, the regularization value is bounded above by a residual term determined by the Lipschitz constant K, which can be minimized by including spectral normalization in G, the norm of the ground truth signal 122 and the RIP constant of the sensing matrix: limg- of 1)) I K'|X**| [1 + A AM )- This bound quantifies the asymptotic regularization performance of the NPN prior in terms of the approximation quality of the learned model and the sensing operator geometry.
Our theoretical analysis has focused on integrating the proposed NPN regularizer into projected- gradient-descent-bale PnP methods. However, the same regularization strategy can be applied to a broad class of learning- baaming-based solvers. For example, in pnporithmmunrollling architectures [39], deep equilibrium models [18], and
frameworks [8] one can insert the NPN proximal step alongside the usual network updates, with end-to-end training of all learnable parameters [44, 9]. Although we demonstrate the empirical benefits of NPN within an unrolled network in Section 5, a full theoretical treatment of these extensions is left for future work. We also incorporate the NPN regularization in Deep Image Prior (DIP) [54] framework (see Appendix A.7 for more details) and in two diffusion-based solvers [15], diffusion posterior sampling (DPS) [12], and DifiPIR [65] (see Appendix A.8 for more details).
5 Experiments
The proposed NPN regularization was evaluated in five imaging inverse problems: compressed sensing, super-resolution, computed tomography, single coil MRI, and deblurring. The method was implemented using the PyTorch framework. For the recovery algorithm, we adopt the FISTA solver [1], with deep denoiser [25]. regularization by denoising (RED) [41], and sparsity prior [1], see Appendix A.5. All simulations were performed on an NVIDIA RTX 4090 GPU, with code in
Compressed Sensing: The single-pixe camera (SPC) is used along with the CIFAR- 10 dataset [28]. with 50, 000 images for training and 10, 000 for testing. All images were resized to 32 x 32. The Adam [26] optimizer was used with a learning rate of 5 x 10 H is random binary sensing matrix with m/n = 0.1. is initialized by QR decomposition, with Algorithm in Appendix A.2. Then, G (for which we used ConvNeXt [33]) and S are optimized following Eq. (3).
MRI: We employed the fastMRI knee single-coil MRI dataset [27], which consists of 900 training images and 73 test images of knee MRIs of 320 x 320. The training set was split into 810 images for training and 90 for validation, and all images were resized to 256 256. For G we used U-Net based architecture [42] and trained it for 60 epochs with a learning rate of x 10-4, using the AdamW optimizer [35] with a weight decay of x 10- and batch size of 4. Performance was evaluated at acceleration factors (AF 21 of 4, 8, and 16, using Radial undersampling mask [32].
Deblurring and Super- Resolution: For the deblurring inverse problem, used 2-D Gaussian kernel with a variance a. For these experiments, we used the CelebA [34] dataset resized to 128 x 128, using 8000 images for training and 2000 for testing. Here, employed a U-Net architecture for the network G. In this case, we used the Adam optimizer with a learning rate of x 10-3 and a batch size of 32. The results for the super-resolution case are shown in Appendix A.
Computed Tomography: We evaluated the proposed regularization on the limited-angle CT inverse problem, using 60 out of 180 total views (spaced every 19) under parallel-bean geometry. We train the DM for 1000 epochs with batch size 4 using the Adam optimizer and learning rate 3 x 10 4 We used cosine variance schedule ranging from 81 1 10 to B1 0.02 with T 1000 time steps. The LoDoPaB-CC dataset [30] was resized to 256 x 256 and used for training, in testing, we used 10 test set slices. G is U-Net which was trained for 100 epochs with learning rate 3 x 10-1 and a batch size of 4 using A.
"github. com/yromarioghh/1P
7

---

CIZ with
CTZ with
PaP-NPN PaP-XPN with
PnF with S
PnP with
[N(( with
-N with
1.00
|Sx Sx to
0.:66
10
0.90
Dectttion
Figure 2: PnP-FISTA convergence analysis in CS. (a) Reconstruction error. (b) Null-space prediction error for (red) Initialization S QR()) from Algorithm I, and (blue) Designed with Eq. (3) and m/n - p/n 0.1. In this case, the CIZ from Definition 3 is highlighted light red and light blue. (c) Acceleration ratio of signal convergence; here, the CTZ is defined as the empirical convergence ratio of the proposed solution that is lower than the baseline (black).
Preponed
34
10
PSNK LBD)
20
40
80
100
300
1010
1200
201
Figure 3: Convergence curves for E 5, 10} in deblurring, AF {4, 12} in MRI.
5.1 Convergence Analysis
To confirm numerically the theoretical results obtained in Sec. 4, we show Fig. with 3 general cases: Baseline (black), where only H is used; PnP-NPN with designed S (blue), where G and S are jointly trained following Eq. (3), with A = 0.001, X2 = 0.01; PnP-NPN with QR (red) where the projection matrix $ QR HD) is fixed (for details about QR(-) see Algorithm in Appendix A.2), only trained, and A = X = 0. Fig. 2(a), we plot the error & ||x x x1 showing that both NPN-QR (red) and NPN-designed (blue) decay much faster than the baseline (black) over the first 75 iterations. Fig. 2(c) tracks the projection onto null-space ||Sx* 1x' which steadily decreases where ||N I I ||S(x' x1 whereas the baseline's reprojection error (dashed black) remains high where. The CIZ is defined in Def. 3. Fig. 2(c) shows the per-step acceleration ratio ||x x x1 ==1||//|x* --xx x 111,, where both NPN curves dip well below the baseline. In this scenario, the CIZ is defined as the iterations when the convergence ratio of the proposed method lower than that of the baseline, with the designed achieving the smallest R2 around & Av 50-75, confirming stronger per-iteration eπrо reduction. The results show that empirical improvements in the algorithm convergence are predictable with the CIZ validating Thcorem 1. Figure shows convergence plots in PSNR for MRI with AF {4, 12), and for deblurring and a {2, 5, 10),
where NPN regularization consistently yields higher reconstruction quality and faster convergence denoisers in PnP-ADMM shown in Appendix A.11 for the image deblurring task. For the CS scenario, the selection and the joint optimization in Eq. (3) fundamental, see Appendices A. and A.12 for detailed analysis.
compared to the baseline. Additional results with other state-of-the--
5.2 Visual results
Figure 4 presents reconstruction results for MRI with an AF = 4, CT with acquisition of 30 views of a total of 180, and for deblurring with 2 using PnP with the DnCNN prior [63]. The estimate & is obtained via equation (1), while the estimate & is obtained with NPN regularization through equation (2). Results show that the learned prior effectively approximates the true nullspace G(y) 2 Sx" with an estimation PSNR of 28.1 18.1 dB in deblurring. 39.13 dB for MRI. Moreover, the reconstruction & obtained with NPN regularization preserves high-frequency details and provides overall improved performance than Xo.

---

51M H
BNMM
MRT -) Deluriin
L715
5313
111
SIM
Figure 4: Deblurring and MRI reconstruction results for PnP and PnP-NPN using DrCNN prior, with dB SNR measurement noise. The measurements are denoted by y. the nonlinear approximation of the null-space of H by y.. and the true null-space by y The estimation Xo is obtained with (1), and & is obtained with the proposed regularization with (2). The ground truth signal is x"
Inverse Problem
Baseline
NPN
NSN
Sparsity
PnP
RED
Sparsity
PnP
RED
DNSN (45]
DDN-C
DDN-I
CS ( = 0.1)
15.93
20.04
17.45
16.15
21.12
17.53
20.10
20.03
20.7
= 4)
36.86
35.99
36.00
38.16
38.08
38.07
35.2
33.7
33.2
Deblurring 29.27
30.78
32.84
31.77
31.42
33.67
33.07
33.03
32.70
Table State-of ofttre-rttttttt comparison for CS, MRI, and Deblurring. For each task, the best results are
highlighted in bold teal, while the second-best results are shown in underline orange
5.3 State-of-the-art comparison
We compared our approach with other state-of-the-art recovery methods that exploit the null-space structure. Mainly, we compare with the DNSN [45] and the DDN-I and DDN-C [6] (For more details the implementation of these methods, refer Appendix A.9). The results for CS, MRI, and Deblurring are in Table 1, showing consistent improvement or competitive performance of the proposed method. In Appendix 1.14, we show a comparison where the same neural network as G. instead of estimating Sx it estimates directly and incorporates it as regularizer in PnP-ADMM algorithm, where we show that the proposed method achieves superior performance.
Method
p/n
CIFAR TO
STLIO
5.4 Performance in data-driven models and dataset generalization
PnP
Unrolling
Unrolling
Baseline
0.0
20.04
24.32
20.09
18.35
In Table we report PSNR (dB) for both PnP and unrolling solvers CIFAR-10 [28] (in-distribution)
0.1
21.12
28.53
19.91
19.64
0.3
21.07
28:75
21.14
20.23
0.5
20.78
18.76
and STL10 [13] (out-of-distribution across projec- tion ratios p/n. n. PnP achieves dB boost at p/n = 0.1 and sustains dB gains on STL10, peak- ing when p/n 0.3. Unrolling delivers up to 5.6 in in-distribution improvement but only 1.9 dB cross- dataset gain at p n 0.3 before declining. Overall, PnP ensures stable generalization, while unrolling maximizes peak PSNR; p/n Rs 0.3 provides the best balance between accuracy and robustness. In general, NPN improves both PnP and unrolling reconstruc- tion performance regarding dataset changes.
0.7
20.09
26.73
20.31
18:45
0.9
20.41
29.90
21.02
19.48
Table 2: Dataset generalization results for SPC in PnP and Unrolling, Each S, G* and Unrolling were optimized with CIFAR-10, and tested with the CIFAR-I0 and STL10. For each dataset, the best results are high- lighted in bold teal while the second-best results are shown underline orange
5.5 Deep Image Prior
We consider the deblurring task with a kernel bandwidth of a 4.0. We train G with the Places365 dataset [64], where we used 28.000 images for training and 7000 images for testing. All images were resized to 128 x 128. The network K was trained following (15) using the Adam optimizer with learning rate of le for 1000 iterations. The network K is U-Net of the same size as G. In Fig. shown the reconstruction performance of DIP and NPN-DIP for different values 7. The results show significant improvements of up to 5 dB and convergence improvements.

---

5.6 Diffusion model solvers
7
NPN-DPS
NPN-DiffPIR
We integrated the proposed regularization term into two widely adopted DM frameworks, DPS [12] and DiftPIR [65]. Table shows the obtained results for
0.0 (Base)
28.22
31.30
10-5
28.55
31.88
different values of For DPS, the NPN regularization consistently improves reconstruction performance by up to 1.85 dB. For DifiPIR, it yields improvements of up to 0.61 dB. Figure shows reconstruction results using NPN within DPS and DifiPIR; yellow arrows highlight improved image details with NPN. Details on the implementation of NPN into DPS and DiAPII are provided in the Appendix A.8.
10-4
28.30
31.91
0.1
30.06
28.98
0.2
30.07
28.57
Table 3: Ablation over y for two meth- ods. Best results are bold teal second- best are underline orange
NPN-DPS
22.5
SSIM PSNR
0.784 [29
0,818 31.72
20.0
17.5
PSNR (R8)
15.0
12.5
Buse DIP
10.0
0.001
DHIPIR 0.872 33.277
NPN-DiLPIR 0.866 33.67
0.01
01
5.0
7 - 0.5
7=
200
600
1000
Boratioo
Figure 5: Performance of DIP and NPN-DIP for different values y in Deblurring.
Figure 6: Reconstruction results using NPN within DPS and DitIPIR.
Limitations
While the proposed method introduces some additional complexity, such as training dedicated neural network for each sensing configuration (H, the networks are lightweight and tailored to specific inverse problems, keeping computational demands modest (see Appendices A.5 and A. 1.14 for details). The current integration into learning-bass based reconstruction frameworks (e.g., unrolling networks) involves two-stage training process, first for the NPN regularizer, then for the reconstruction model, but this modular design enables flexible adaptation. Future work could explore joint end-to-end training of the NPN regularizer and reconstruction network to further enhance performance and efficiency. Unlike NSN-based methods [45, 6), which directly reconstruct, our approach learns subspace projection offline, adding an extra training step but offering improved compatibility with different solvers. Although we devise some design criteria for selection S, there are scenarios in which the method can fail upon this design, we provide detailed discussion in A. on this aspect from the point of view of finding non-linear relations between Hx and Sx.
7 Conclusion and future outlooks
We introduce Non-Linear Projections of the Null-Space for regularizing imaging inverse problems. Intuitively, the regularization promotes selective coefficients of the signal in the null-space of the sensing matrix. This formulation allows ficxibility in what features the null-space we can exploit. Our proposed method is equipped with strong theoretical guarantees for plug-and-play restoration algorithms, showing that the proposed regularization has a zone of convergence improvement con- trolled by the network error. Additionally, we show that our regularizer converges to constant depending also on the network estimation error when the algorithm reaches the optimum. We validate theoretical findings in five distinctive imaging inverse problems: compressed sensing, magnetic resonance, super-resolution, computed tomography, and deblurring. Results validate the theoretical developments, and we have improved performance with state-of-the-an methods. This approach opens new frontiers to regularize the imaging inverse problems in different solvers, such as deep equilibrium models and consensus equilibrium formulations.
10

---

Acknowledgements
This work was supported in part by the Agencia Nacional de Hidrocarburos (ANH) and the Ministerio de Ciencia, Tecnologia Innovación (MINCIENCIAS), under contract 045-2025, and in part by the Army Research Office/Laboratory under grant number W911NF-25-10015, VIE from UIS project 8087. The views and conclusions contained in this document are those of the authors and should not be interpreted as representing the official policies, either expressed or implied, of the Army Research Laboratory or the U.S. Government.
References
[1] Amir Beck and Mare Teboulle. A fast iterative shrinkagg-threshoolinng algorithm for linear inverse problems. SIAM journal on imaging sciences, 2(1):183-202, 2009.
[2] Antoni Buades, Bartomeu Coll, and Jean-Michel Morel. Non-Local Means Denoising. Image Processing On Line,
2011.
[3] Harold C Burger, Christian J Schuler, and Stefan Harmeling. Image denoising: Can plain neural networks compete with bm3d? In 2012 IEEE conference computer vision and pattern recognition, pages 2392-2399, IEEE, 2012.
[4] E. J. Candes and M. B. Wakin. An introduction to compressive sampling. Processing Magazine, 25(2):21-30, 2008.
IEFE Signal
[5] Bin Chen, Jiechong Song, Jingfen Xie, and Jian Zhang. Deep physics-guided unrolling generalization for compressed sensing. International Journal of Computer Vision, :1:(11) 2887, 2023.
[6] Dongdong Chen and Mike E Davies. Deep decomposition learning for inverse imaging problems. In Computer Vision-ECCV ECCV 2020: 16th European Conference, Glasgow, UK, August 23-28, 2020, Proceedings, Part XXVIII 16, pages 510-526. Springer, 2020.
[7] Dongdong Chen, Julián Tachella, and Mike E Davies. Equivariant imaging: Learning beyond the range space. In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 4379 4388, 2021.
[8] Tianlong Chen, Xiaohan Chen, Wuyang Chen, Howard Heaton, Jialin Liu, Zhangyang Wang, and Wotao Yin. Learning to optimize: A primer and a benchmark. Journal of Machine Learning Research, 23(189): 2022.
[9] Xiaohan Chen, Jialin Liu, Zhangyang Wang, and Wotao Yin. Theoretical linear convergence of unfolded ista and its practical weights and thresholds. Advances in Neural Information Processing Systems, 31, 2018
[10] Yurong Chen, Yaonan Wang, and Hui Zhang. Unsupervised range-nullspace learning prior for multispectral images reconstruction. IEEE Transactions on Image Processing, 2025.
[11] Xinhua Cheng, Nan Zhang, Jiwen Yu, Yinhuai Wang, Ge Li, and Jian Zhang. Null-space diffusion sampling for zero-shot point cloud completion. In IJCAL, pages 618 0--66, 2023.
[12] Hyungjin Chung, Jeongsol Kim, Michael T Mccann, Mare L Klasky, and Jong Chul Ye. Diffu- sion posterior sampling for general noisy inverse problems. arXiv preprint arXiv:2209.14687 2022.
[13] Adam Coates, Andrew Ng. and Honglak Lee. An analysis of single-layer networks in unsuper- vised feature learning. In Proceedings of the fourteenth international conference on artificial intelligence and statistics, pages 215-222 JMLR Workshop and Conference Proceedings, 2011.
[14] Kostadin Dabov, Alessandro Foi, Vladimir Katkovnik, and Karen Egiazarian. Image denoising with matching and 3d filtering. In Image processing: algorithms and systems, neural networks, and machine learning, volume 6064, pages 354-365 SPIE, 2006.
11

---

[15] Giannis Daras, Hyungjin Chung, Chich-Hsin Lai, Yuki Mitsufuji, Jong Chul Ye, Peyman Milanfar, Alexandros G Dimakis, and Mauricio Delbracio. A survey on diffusion models for inverse problems. arXiv preprint arXiv:2410.00083, 2024.
[16] Eduardo S. L. Gastal and Manuel M. Oliveira. Domain transform for edge-aware image and video processing. In ACM SIGGRAPH 2011 Papers, SIGGRAPH '11, New York, NY, USA, 2011. Association for Computing Machinery.
[17] Tatiana Gelvez, Hoover Rueda, and Henry Arguello. Joint sparse and low rank recovery algorithm for compressive hyperspectral imaging. Appl. Opt., 56(24))6785 Aug 2017.
[18] Davis Gilton, Gregory Ongie, and Rebecca Willett. Deep equilibrium architectures for inverse problems in imaging. IEEE Transactions Computational Imaging, 7:1123-118 2021.
[19] Gene H Golub, Per Christian Hansen, and Dianne P O'LLery. Tikhonov regularization and total least squares, SIAM journal on matrix analysis and applications, 21(1):185-194, 1999.
[20] Romario Homudróón-Hurtdde Roman Jacome, Leon Suarez, Emmanuel Martinez, and Henry Arguello. Improving compressive imaging recovery via measurement augmentation. In ICASSP 2025 2025 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). pages 1-5, 2025.
[21] Bahadir Gunturk and Xin Li. Image restoration. CRC Press, 2018.
[22] Jonathan Ho, Ajay Jain, and Pieter Abbeel. Denoising diffusion probabilistic models. In
H. Larochelle, M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin, editors, Advances in Neural Information Processing Systems, volume 33, pages 6840-6851. Curran Associates, Inc., 2020.
[23] Yuyang Hu, Mauricio Delbracio, Peyman Milanfar, and Ulugbek Kamilov. A restoration network as implicit prior. arXiv preprint 2rXiv:2310.01391, 2023.
[24] Yuyang Hu, Albert Peng, Weijie Gan, Peyman Milanfar, Mauricio Delbracio, and Ulugbek Kamilov. Stochastic deep restoration priors for imaging inverse problems. arXiv preprint arXiv:2410.02057, 2024.
[25] Kamilov et al. Plug-and-play methods for integrating physical and learned models in compu- tational imaging: Theory, algorithms, and applications. IEEE Sig. Proc. Mag.. 40(1):85-99) 97, 2023.
[26] Diederík P Kingma and Jimmy Ba. Adam: A method for stochastic optimization. arXiv preprint 4rXiv:1412.6980, 2014.
[27] Florian Knoll, Jure Zbontar, Anuroop Sriram, Matthew J Muckley, Mary Bruno, Aaron Defazio, Mare Parente, Krzysztof J Geras, Joe Katsnelson, Hersh Chandarana, Zizhao Zhang, Michal Drozdzalv, Adriana Romero, Michael Rabbat, Pascal Vincent, James Pinkerton, Duo Wang, Nafissa Yakubova, Erich Owens, C Lawrence Zitnick, Michael P Recht, Daniel Sodickson, and Yvonne W Lui. FastMRI: A publicly available raw k-space and DICOM dataset of knee images for accelerated MR image reconstruction using machine learning. Radiol. Artif. Intell., ((1)ee19000) January 2020.
[28] Alex Krizhevsky, Geoffrey Hinton, et al. Learning multiple layers of features from tiny images. 2009.
[29] Joseph Kuo, Jason Granstedt, Umberto Villa, and Mark A Anastasio. Learning projection operator onto the null space of linear imaging operator. In Medical Imaging 2021: Physics of Medical Imaging, volume 11595, pages 1019-1025, SPIE, 2021.
[30] Johannes Leuschner, Maximilian Schmidt, Daniel Otero Baguer, and Peter Maass. LoDoPaB- CT, benchmark dataset for low-dose computed tomography reconstruction. Sei. Data, 8(1):109, April 2021.
[31] Jiaming Liu, Salman Asif, Brendt Wohlberg, and Ulugbek Kamilov. Recovery analysis for plug-and-play priors using the restricted eigenvalue condition. Advances Neural Information Processing Systems, 34:5921-5933, 2021.
12

---

[32] Jing Liu and David Saloner. Accelerated MRI with CIRcular cartesian UnderSampling (CIR- CUS): variable density cartesian sampling strategy for compressed sensing and parallel imaging. Quant. Imaging Med. Surg., 4(1):57-67, February 2014.
[33] Zhuang Liu, Hanzi Mao, Chao- Wu, Christoph Feichtenhofer. Trevor Darrell, and Saining Xie. A convnet for the 2020s. In Proceedings of the IEEE/CVF conference on computer vision and pattern recognition, pages 1 1196--19986 2022.
[34] Ziwei Liu, Ping Luo, Xiaogang Wang, and Xiaoou Tang. Deep learning face attributes in the wild. In Proceedings of International Conference on Computer Vision (ICCV), December 2015.
[35]
Iyy Loshchilov and Frank Hutter. Decoupled weight decay regularization, 2019.
[36] Michael Lustig, David L. Donoho, Juan M Santos, and John M Pauly. Compressed sensing mri.
IEEE signal processing magazine, 25(2):72-82, 2008.
[37] Lei Ma, Hulin Kuang, Jin Liu, Chengchao Shen, and Jianxin Wang. Improving medical image denoising via a lightweight plug-and-playy module. In 2023 IEEE International Conference on Bioinformatics and Biomedicine (BIBM), pages 1350-1355. IEEE, 2023.
[38] Emmanuel Martinez, Roman Jacome, Romario Gualdrón Hurtado, Iñaki Esnaola, and Henry Arguello. Compressive sensing with augmented measurements via generative self-distillation. In 2025 IEEE Statistical Signal Processing Workshop (SSP), pages 31-35, 2025.
[39] Vishal Monga, Yuelong Li, and Yonina C Eldar. Algorithm unrolling: Interpretable, efficient deep learning for signal and image processing. IEEE Signal Processing Magazine, 38(2): 2021.
[40] Neal Parikh and Stephen Boyd. Proximal algorithms. Foundations and Trends in optimization, 13)::127-23, 2014.
[41] Yaniv Romano, Michael Elad, and Peyman Milanfar: The little engine that could: Regularization 2017.
by denoising (red). STAM Journal Imaging Sciences,
[42] Olaf Ronneberger, Philipp Fischer, and Thomas Brox. U-net: Convolutional networks for biomedical image segmentation. In Medical Image Computing and Computer-Assister Intervention-Men MICCAI 2015: 18th International Conference, Munich, Germany. October 5-9. 2015, Proceedings, Part III 18, pages 234-241. Springer, 2015.
[43] Ernest Ryu, Jialin Liu, Sicheng Wang, Xiaohan Chen, Zhangyang Wang, and Wotao Yin. Plug-and-play methods provably converge with properly trained denoisers. International Conference on Machine Learning, pages 5546-5557. PMLR, 2019.
[44] Jonathan Scarlett, Reinhard Heckel, Miguel RD Rodrigues, Paul Hand, and Yonina C Eldar. Theoretical perspectives on deep learning methods in inverse problems. IEEE journal on selected areas in information theory, 3(3):433 453, 2023.
[45] Johannes Schwab, Stephan Antholzer, and Markus Haltmeier. Deep null space learning for inverse problems: convergence analysis and rates. Inverse Problems, 35(2):025008, 2019.
[46] Parikshit Shah and Venkat Chandrasekaran. Iterative projections for signal identification on manifolds: Global recovery guarantees. In 2011 49th Annual Allerton Conference on Communication, Control, and Computing (Allerton), pages 760-767. IEEE, 2011.
[47] Jascha pervised learning using nonequilibrium thermodynamics. In Francis Bach and David Blei. editors, Proceedings of the 32nd International Conference on Machine Learning, volume 37 of Proceedings of Machine Learning Research, pages 2266-2265, Lille, France, 07-09 Jul 2015.
Eric Weiss, Niru Maheswaranathan, and Surya Ganguli. Deep unsu-
PMLR.
[48] Bowen Song. Soo Min Kwon, Zecheng Zhang, Xinyu Hu, Qing Qu, and Liyue Shen. Solv- ing inverse problems with latent diffusion models via hard data consistency. arXiv preprint arXiv:2307.00 08 123, 2023.
13

---

[49] Jiaming Song. Chenlin Meng, and Stefano Ermon. Denoising diffusion implicit models. In International Conference on Learning Representations, 2021.
[50] Yang Song, Jascha Sohl-Dickstein, Diederik P Kingma, Abhishek Kumar, Stefano Ermon, and Ben Poole. Score-based generative modeling through stochastic differential equations. In International Conference on Learning Representations, 2021.
[51] Afonso M Teodoro, José M Bioucas-Dias, and Mário AT Figueiredo. Image restoration and reconstruction using targeted alug-and-play priors. IEEE Transactions Computational Imaging, 5(4):675-686, 2019.
[52] Matthieu Terris, Ulugbek S Kamilov, and Thomas Moreau. Fire: Fixed-points of restoration priors for solving inverse problems. arXiv preprint arXiv:2411.18977 2024.
[53] Matthieu Terris, Thomas Moreau, Nelly Pustelnik, and Julian Tachella. Equivariant plug-and- play image reconstruction. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 25255 252664, 2024.
[54] Dmitry Ulyanov, Andrea Vedaldi, and Victor Lempitsky. Deep image prior. In Proceedings of the IEEE conference computer vision and pattern recognition, pages 9446-9454 2018.
[55] Singanallur V. Venkatakrishnan, Charles A. Bouman, and Brendt Wohlberg. Plug-and-play pri- ors for model based reconstruction. In 2013 IEEE Global Conference on Signal and Information Processing, pages 945-948, 2013.
[56] Yinhuai Wang. Yujie Hu, Jiwen Yu, and Jian Zhang. Gan prior based null spll-space learning for consistent super-resolution. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 37, pages 2724 2732, 2023.
[57] Yinhuai Wang, Jiwen Yu, Runyi Yu, and Jian Zhang. Unlimited-size diffussonn restoration. In Proceedings of the IEEE/CVF conference computer vision and pattern recognition, pages 1160-1167, 2023.
[58] Yinhuai Wang, Jiwen Yu, and Jian Zhang. Zero-shot image restoration using denoising diffusion null-space model. arXiv preprint arXiv:2212.00490 2022.
[59] Martin J Willemink, Pim A de Jong, Tim Leiner, Linda M de Heer, Rutger AJ Nievelstein, Ricardo PJ Budde, and Arnold MR Schilham. Iterative reconstruction techniques for computed tomography part 1: technical principles. European radiology, 163,,,1163,, 2013.
[60] Zamir et al. Restormer: Efficient transformer for high-resolution image restoration. In IEFF/CVF CVPR, pages 5728-5739, 2022.
[61] Zhiyuan Zha, Bihan Wen, Xin Yuan, Saiprasad Ravishankar, Jiantao Zhou, and Ce Zhu. Learning nonlocal sparse and low-rank models for image compressive sensing: Nonlocal sparse and low-rank modeling. IEEE Signal Processing Magazine, 40(1):32-44 44, 2023.
[62] Kai Zhang, Yawei Li, Wangmeng Zuo, Lei Zhang, Lue Van Gool, and Radu Timofte. Plug-and- play image restoration with deep denoiser prior. IEEE Transactions on Pattern Analysis and Machine Intelligence, 44(10):6360-6376 2021.
[63] Kai Zhang, Wangmeng Zuo, Yunjin Chen, Deyu Meng, and Lei Zhang. Beyond gaussian denoiser: Residual learning of deep enn for image denoising. IEEE transactions on image processing, 26(7):3142 3155, 2017.
[64] Bolei Zhou, Agata Lapedriza, Aditya Khosla, Aude Oliva, and Antonio Torralba. Places:
10 million image database for scene recognition. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2017.
[65] Yuanzhi Zhu, Kai Zhang, Jingyun Liang, Jiczhang Cao, Bihan Wen, Radu Timofte, and Luc Van Gool. Denoising diffusion models for alug-and-play image restoration. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR) Workshops, pages June 2023.
14

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: The abstract and introduction clearly states the technical contributions of the paper along with the description of experimental validation and the relation with prior art.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
It is fine to include aspirational goals as motivation long is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: We included a dedicated section on limitations to thoroughly discuss the constraints and potential weaknesses of our approach.
Guidelines:
The answer NA means that the paper has limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create a separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (c.g.. independence assumptions, noiseless settings, model sell-speciiiatiinn asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g. if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect on the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-tex, system might not be used reliably to provide closed captions for online lectures because it fails handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and faimesss.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community, Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
15

---

Answer: [Yes]
Justification: All theoretical results in the paper are accompanied by clearly stated assump- tions and complete proofs. Each theorem, definition, and assumption is properly numbered, cross-referenced, and supported by formal derivations, either the main paper or in the supplementary material, with intuitive sketches provided where appropriate.
Guidelines:
The answer NA means that the paper does not include theoretical results. All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced in the statement of any theorems.
The proof' can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should properly referenced.
4. Experimental result reproducibility
Does the paper fully disclose all the information needed to reproduce the main ex- perimental results the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Question:
Answer: [Yes]
Justification: The paper provides a comprehensive experimental description with its corre- sponding settings for reproducing the results
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurlPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make it clear how to reproduce that algorithm.
(b) the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (e.g., a large language model), then there should either be way access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (c.g., to registered users), but it should be possible for other researchers to have some path to reproducing verifying the results.
16

---

5. Open access to data and code
Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
Question:
Answer: [Yes]
Justification: The code will be uploaded in the supplementary material with sufficient instructions to reproduce the results of the paper.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nips. cc/
Mabllic/guides//odSSttaaae sot.esuubbbssiioooooooooy. for more details.
While we encourage the release of code and data, we understand that this might not be possible, "No" an acceptable answer. Papers cannot be rejected simply for including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https: / //niss. cc/ccc/ppu wwwpcppwpcc//eeeecceeeceeeeeetle gttcc//gg/////eeeeeeeecc
for more details.
The authors should provide instructions on data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental setting/details
Does the paper specify all the training and test details (c.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Question:
Answer: [Yes]
Justification: The paper provides comprehensive details regarding the experimental setup, allowing for reproducibility and clear understanding of the reported results.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core the paper to a level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, as supplemental material.
7. Experiment statistical significance
Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Answer: [No]
Justification:
The paper reports only the mean performance metrics for all experiments.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by emror bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
17

---

The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
It should be clear whether the error bar is the standard deviation the standard error of the mean.
It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful to show in tables or figures symmetric error bars that would yield results that are out of range (c.g. negative error rates).
If error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of exccution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification: The paper mentions the computational resources employed obtain the experimental results, and it complemented with the memory and time of execution the supplementary material
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported the paper (e.g., preliminary or failed experiments that didn't make into the paper).
9. Code of ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https: : http//eeurip...ee cc/publicc/thiccauiieeee
Answer: [Yes]
Justification: This research fully complies with the NeurIPS Code of Ethics. does not involve any sensitive information, personal data, human participants, applications with significant societal risks. All datasets utilized are publicly accessible and ethically obtained, and every experiment was conducted transparently, fairly, and reproducibly.
Guidelines:
The answer NA means that the authors have not reviewed the NeurlPS Code of Ethics.
If the authors answer No, they should explain the special cireumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., if there is special consid- cration due to laws or regulations their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative socictal impacts of the work performed?
Answer: [Yes]
18

---

The paper considers applications in medical imaging (MRI) which opens new frontiers for further automated and precise medical diagnosis as well as improved image restoration task for commercial cameras.
Justification:
Guidelines:
The answer NA means that there no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation. generating fake profiles, surveillance), faimess considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology.
If there negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data or models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: The generated weights, models, and results do not have risk for misuse as only a public dataset was used.
Guidelines:
The answer NA means that the paper poses such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards challenging, and many papers do not require this, but encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
The paper properly cites original work on the used models, data and technical developments.
Justification:
Guidelines:
The answer NA means that the paper does not use existing assets.
19

---

The authors should cite the original paper that produced the code package or dataset.
The authors should state which version of the asset is used and, if possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset. For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswi thcode. com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. If this information is available online, the authors are encouraged to reach out to the asset's ss creators.
13. New assets
Question: Are new assets introduced in the paper well documented and the documentation provided alongside the assets?
Answer: [Yes]
Justification: We release both code and trained models as part of this work. All assets are thoroughly documented, including setup instructions, usage examples, model architecture details, training procedures, and licensing information. The assets are provided in the anonymized repository for review support reproducibility and ease of use.
Guidelines:
The answer NA means that the paper does not release assets.
Researchers should communicate the details the dataset/code/model as part their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset used.
At submission time, remember anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification:
The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According to the NeurlPS Code of Ethies, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals or equivalent for rescarch with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
20

---

Answer:
[NA]
The paper does not involve crowdsourcing nor research with human subjects Guidelines:
Justification:
The answer NA means that the paper does involve crowdsourcing nor research with
human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent)
may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurlPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such the institution conducting the review.
16. Declaration of LLM usage
Question: Does the paper describe the usage of LLMs if it is an important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, or formatting purposes and does not impact the core methodology, scientific rigorousness, originality of the research, declaration is not required.
Answer: [NA]
Justification: This research does not involve LLMs as any important, original, or non- standard components.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, non non- non-standard components.
Please refer to our LLM policy (https for what should or should not be described.
21

---

A Technical Appendices
A.1 Settings for experiment in Figure
To illustrate this concept, we developed a toy example with x" R. H. E g2x3, and constructed an orthonormal vector S 1 relative to the rows of H, following Algorithm 1. A training dataset p(x) was generated, consisting of random points within circle centered at the origin. Using this dataset, we trained two-layer neural network defined as G(Hx) = Vo WHx), with parameters W E Rxxx2 and V E 1 with k = 50 hidden neurons, optimized via the problem stated in (3). For comparison, we also trained a network of similar size to directly reconstruct Xo from measurements Hx".
A.2 Algorithms for designing
We developed algorithm for obtaining S from H. satisfying that SH -. The algorithm is based on the QR decomposition, first, computing a full QR decomposition of H T E =R=*",, yielding an orthonormal basis Quiil Raxn for IR, The columns from 1 to of Qruil form basis for (H), which we denote by N E Bax(n-m)) To construct a subspace of the null space, the algorithm samples a random Gaussian matrix P E (m=m))xp, which is orthonormalized via QR decomposition to produce U E ((- Xp, This ensures that the resulting subspace is both diverse and well-conditioned. Finally, the matrix S is obtained
N E R2*n which consists of p orthonormal vectors that span a random p-dimensional subspace within Nu)l(H).
Algorithm 1 GENERATE ORTHONORMAL Rows TO H VIA QR DECOMPOSITION
Require: Matrix H E Bm desired number of rows
Ensure: Matrix S E TR=x" whose rows orthonormal and lie in Nu)l(H)
I: Quall t- QR(HT)
2: N Quai(:
1:m]
Nullspace basis, size (n m)
3: Sample P ~ 1) (nnm))Xp
4: + QR(P)
((-mm) with orthonormal columns Resulting x matrix of orthonormal rows
5: S U N
return S
A.3 Proof of Theorem 1
Here we provide the technical proof for Theorem 1.
Theorem (PnP-NPN Convergence). Consider the fidelity g(x)= = ||H - ylž the denoiser Drr satisfies assumption 2, the matrix S satisfies the RIP condition 2 with constant Ai,, E [0, 1). The optimized network G satisfies assumption 3. We run iterations as in (4) in the convergence improvement zone for e == 1...., Lg. The algorithm converges to a fixed point if
(1 + (|| a(H H SST + + + A I1 <
(9)
Proof:
||x 4+
111 = -
Assmpption 2
(1+ ||x a (H (Hx' yS - (Sx" + +N( ())x)))) - -12
= (1+ (0+0))((1
||[[ aI H S
(xx
triang, ineg.
(1+ 1- a H 8) | - (1 8) S N(Hx*)]
Defuition I
(1+0) || I H
s)||
| -x-1| + (1 S||SS (x' - -x.) ||
Deficition
(1+ (
H H />
(1+ - ) 1| - x
22

---

A.4 Proof of Theorem 2
Theorem 2 (Convergence of NPN Regularization), Consider the iterations of the PnP-NPN algorithm defined in (4) for 8 = 1. L. Assume that the estimation error of the trained network G* satisfies Assumption 1, and that both S and H satisfy the Restricted Isometry Property (RIP) over the manifold MD with constants AMp: AMD E (0, 1). Further assume that the residual term N(Hxx) ) is K-
Lipschitz continuous, thus [[N(Hx")II < K(1 + A Mp }|xx")] Then, after E iteeationn, the NPN
regularization term satisfies the following upper bound:
|[G(y) - -SX++1| - C1|lx
x'11 K||x* (1(1 + ADD) + C2lxx - xl
(10)
where
C1 (1+
2 (1 (1+ + K(1 + A Mp )1 + + C2 =
MDD
1
/2a
Proof: We begin by recalling the definition of the regularization function:
Φ(x) |Sx Sx' N(Hx*)II
From the iteration difference, we have:
((x++)= - IS(x*** x3]12 2(x x S (ST(Sx -Sx" -X""-X(Hx))) (II)
Define the intermediate step:
as q = - (Sxx - Sxx NN(Hx)),, 0.0.
Optimality implies:
|x - q'll S ||x 9'1l
Substituting the definition of q' and rearranging gives:
||x[+1 x" + a T (Sx² Sx" - N(Hx < |xx - x + as T - Sx" - -N(xx*)) ?||?
Expanding both sides and reorganizing terms, we obtain:
20(x*t1 ++1 xx, S - Sx" - N(Hxx))) I ||x x112 1|x - x112
+ 2a(xx - x', (Sx -Sxx - N(Hx))))
Dividing by 2a and substituting back into (11) yields:
o(x x(xx) < ||S(x ++1 xf )11
+ ||x x1122
1
1
||x e+1 - x x11
20
2a
-x" ST SSx' Sxx 1**-NNIU -N(Hx")))
(12)
Next, applying Assumptions 2 and 3, we bound:
(x' x TSx' Sxx - N(Hx*)))
< |S(x xx'|| I12 + K(1 + A (1 + A
Mp
||x*1|] |11x x"1
Substituting this back into (12), we get:
d(x(t)) {xx)) S ||(xx+11 -x)011
1
1x
e+ - x²11
2v
+ 2a
1x x112 1|5(x xx1)112
+ K(1 + Mp (1 + A )||x" 11||= x"11
(13)
Using the triangle inequality and the RIP property again, we simplify:
o(x²)= S (1 + A )x
x'1| +|]++K[1 + AMD ||x
23

---

Hence, regrouping terms, we achieve final compact form analogous to the desired inequality:
[[G(y) Sx 0(x(+1)
$ 2a
(1+ ) K(1 A M (1+
AMD )1xe))) ||x
x'1
K(1 ( || + (1+ /20
1
-
x
(14)
A.5
Neural network details
For the single-pixel camera experiments, we employ ConsNeXt-inspired backbone [33]. The network has five ConvNeXt blocks, each comprising two successive
SPC:
convolutions with ReLU activations. A final 3 3 convolutional layer projects the 128 features back to one channel. The output feature map is flattened and passed through linear module to match the measurement dimensionality. We use cosine-based positional encoding with dimension of 256, and set the number of blocks to 5 and the base channel width to 128.
MRI and Deblurring: To train G for the MRI and Deblurring experiments, we use U-Net ar- chitecture [42] with three downscaling and three upscaling modules. Each module consists of two consecutive Conv + ReLU blocks. The downscaling path uses max pooling for spatial re- duction, starting with 128 filters and increasing up to 1,024 the bottleneck. The upscaling path performs nearest-neighhor interpolation before each module, progressively reducing the number of filters. final 2D convolutional layer without activation produces the output. Skip connections link corresponding layers in the encoder and decoder.
A.6 Plug and Play algorithms
This work uses the PnP-FISTA, its unrolled version, and RED-FISTA algorithms to validate the proposed approach. Below is shown algorithms for PnP and RED formulation and their NPN- regularized counterparts. For unrolling FISTA, the only change is that the denoiser D, changed to trainable deep neural network that is optimized in an end-to-end manner.
Algorithm 2 PnP-FISTA
Algorithm 3 NPN-PAP-FISTA
Require: L, H I: x" = z 0, f
Require: L, H.H.S,yCC. ,
1: x z = 0, t=1
2: for = L,...: L do
2: for l
L do
x - aH (He' y)
3: x' + z' HH (Hz' y) + 7 -Sx1-G'(x))
4: x" + (x*)
4:
x" (x
5:
I'
== It 1++1(1
1= ++11+
2
7: 2t11 t x² + c' (x* x(-1)
7:
z 241
x" + ('- (x'-
8: end for
end for return
9: return '
Algorithm RED-FISTA
Algorithm 5 NPN-RED-FISTA
Require: L, H. a 1: x" z' = 0, t
Require: L. H.H.S,y,C ,a,,,
a,,
1: x" z = 0, t=1
2: for l,...,K do x + -
2: for l =
T. do
aH
x" + - (H" (Hz' (H'H'
- G"
(x))))
Mlz (x'r))
4: =t=
1'
1= /14111
l = i+
5:
6: 2t1 x + d-1 (xf -1)
6:
z [+1 + x + (' (x'xx-xx
7: end for
7: end for
8: return '
return x
24

---

A.7 Deep Image Prior
Deep Image Prior (DIP) [54] is an unsupervised reconstruction framework that leverages untrained neural networks to reconstruct the underlying signal. This approach, since it only compares with the measurements y. its solution may overfit to the noise of the measurements. Therefore, the use of the NPN provides a suitable approach to improve image consistency reconstructions. The optimization problem with the proposed NPN regularization is
K* arg min y
+ Sk(z)l
==) z~ ~~,(011), (15)
where y is a regularization parameter.
A.8 Diffusion-based Solvers
Diffusion models (DMs) [47, 22, 49, 50] have recently gained attention due to their exceptional capability in modeling complex image distributions via an iterative noising-denoising process. Condi- tioning DMs entails guiding their generative reverse diffusion process using measured data to ensure reconstructions align with the measurements. We have integrated our proposed regularization term into two widely adopted diffusion- model frameworks, Diffusion Posterior Sampling (DPS) [12] and DifIPIR [65]. These frameworks serve as canonical diffusion pipelines upon which newer methods build. Our regularization could likewise be incorporated into other approaches, such as latent-space diffusion models 48).
DPS [12]: we denote N is the number of reverse diffusion steps, and i E O..... N - 1 is the reverse- time index; x; E R" the current latent state and x (O I) is the Gaussian start; so (x., is the score/noise estimate produced by the network with parameters 0; Xo is the network's s prediction of the clean sample at step ag E (0, is the per-step retention factor, B
a the noise increment, and a, = I ay is the cumulative product (with ã = 1); Ç. > is the data-consistency step size and o; 2 0 is the sampling noise scale at step i; I), is i.i.d. Gaussian noise; - denotes the -oessssaccons iterate before applying the gradient correction; Vx; | H is the gradient of the quadratic data-fidelity term with respect to x; (through the dependence *o(x:));
Algorithm 6 DPS Sampling
Algorithm NPN-DPS Sampling
Require: N. H. y, { ) (ā. ). step
Require: N, H, y, { and noise scales
(ā. JM
b step sizes
sizes and noise scales
I: XN ~ (0,
1: XN (O.
2: for i == N
do
2:
...., do
3: $ solxi, i)
3: 8 + 50 (x i)
score noisc-prediction net Xu + (xi + (1 )8)
4: *o Varg
(x:
aa)s))
4:
5:
~ I0,
5:
z 7~ N (0, I)
X/-1
1-a4
Xi
6:
x/-1 + - 1) 1-a;
x; + Vo 18,
a; 8o + Ō
y Va--TA
T-o;
20 + ō:z
7:
Xi-1 +- x/-1
C/( Vx Ily H 11 +
7: X
+ X - Si Vxx
x 1
7 Gl("(y)) -Soll)
8: end for
8: end for
9: return Xp
9: return Xo
DifIPRR [65]. On > O denotes the standard deviation of the measurement noise, À 0 is the data-proximal penalty that trades off data fidelity and the denoiser prior inside the subproblem; p. X, of. @ is the peration--ependee
weight used in the proximal objective at step i; (i) is the score-model denoised prediction of the clean sample step i (before enforcing data consistency);
(9) is the solution of the data-pron antaippoximml subproblem at step (i.c., the data-consistent refinement
of x ); ē = os) 1/2
(x, 11/2
Vo,, x is the effective noise estimate implied by (x,, x) )
6; ~ A I) is the fresh Gaussian noise injected at step i; Ç E [0, mixes deterministic and stochastic updates in the reverse diffusion (ç = 0 fully deterministic,
fully stochastic).
25

---

Algorithm 8 DiffPIR Sampling
Algorithm 9 NPN-DiffPIR Sampling
Require: N, H, on, fōr fem1 C. A
för N
Require: N, H, y. ou. for J C. X, 7, G" 1: Precompute P: /o for = 1, .V.....V
I: Precompute Pr + / for = II....N
Xx (0,1)
XN ~ (0,11)
3: for = N, N...., do
3: for i = N...... N...., do
4: $ + i)
4: $ + so( i)
5:
x(i)
f vas (xi + (1 ())
5:
x + ,a,
x + (1 - 5,)8)
6:
((0)
f arg minx lly Hxl13 + pillx
6:
%0 t arg minx lly Hx| + +illx
X 12
X
1 11 + 7 || (y) - -xx1
7:
6t
1
Xi -vacx
7:
tt vi-a;
Xi v x
vT-ai
8:
6, ~
N/,11 I)
8:
6, I(0,1)
9:
xi-1
vāi-1 x60
9:
x;-1
v vōi-1 x/1
v1 - di-1 (VI- ci)
v1 - a (VT viei)
10: end for
10: end for
11: return x
11: return x())
A.9 Implementation details on NSN-based methods
Recall the NSN-based models used in comparison:
1. Deep null space network (DNSN)(45): & Hly + (I -I- H'H)((HH)))
2. Deep decomposition network cascade (DDN-C) [6]: H'y + (H(,(
P., (R,,H,y + (H.,))
3. Deep decomposition network independent (DDN--) [6]: Híy (R, +
Pn (R,, (H'y))
We used the source code of [6] to implement the models R. R, and R.
Network R and R, is lightweight version of the U-Net architecture for image segmenta- tion. It features an encoder-decoder structure with skip connections. The encoder consists of five convolutional blocks, each followed by max- -poling to extract features. The decoder upsamples using transposed convolutions and refines the features with additional convo- lutions, leveraging skip connections to combine low- and high-level features. The final output is produced through 1 1 convolution, reducing the output to the desired number of channels. This compact design makes it efficient for tasks with limited resources.
Network R,: It consists of series of convolutional layers arranged
sequential block structure, where each block performs transformation of the input image through convo- lutional operations, followed by activation functions, such as ReLU, and normalization techniques like batch normalization. The model starts with input image, and through several layers of convolutions, the features of the image are progressively refined. Each layer extracts relevant features, while ReLU activations introduce non-linearity to improve the model's capacity to capture complex patterns. Batch normalization layers are added stabilize training and speed up convergence by normalizing the output of each convolu- tional layer. The final output layer of the network reconstructs the denoised image. The architecture is designed with depth of 17 layers.
All models trained for 100 epochs using the Adam optimizer with learning rate of 1e-3 using a mean- squared-error loss function.
A.10 Super resolution experiments
For this scenario, we use the Places365 dataset We set an SR factor SRF = was performed with bilinear interpolation. The forward operator is modeled as DB E xn
V n me
4; downsampling
26

---

Metric
PnP-FISTA Baseline
PnP-FISTA NPN
Estimation of Sx
PSNR SSIM
22.01
23.75
24.95
0.562
0.692
0.697
Table 4: Quantitative comparison of PnP-FISTA and PnP-FISTA NPN with the estimation of Sx.
where B IRux structured Toeplitz matrix implementing Gaussian blur with bandwidth o 2.0
and D E Rmxx is the downsampling matrix.
To recover the information removed by H. define the null-space projection operator, similar to the deblurring case, S i + j] 1 - h[j], i == 1..... m. = 1,..., n, which captures the high-frequency details of the signal. We use the PnP-FISTA algorithm for evaluation with 60 iterations using = 0.5 and y = 0.1. In the following table, we show the reconstruction performance with baseline PnP-FISTA and NPN PnP-FISTA, showing an improvement of 1.74 dB in PSNR and 0.13 in SSIM. Additionally, we show the network estimation of the null-space Sx metrics, showing good estimation of the high-frequency details from the low-resolution image.
A.11 Additional Experiments
In Tables and 6, we provide additional results of MRT and deblurring, respectively, using PnP, RED, and sparsity priors for acceleration factors (AFs) of 4, 12 (with dB SNR measurement noise) and a of 2, 6, 10. These results validate the robustness NPN regularization under different numbers of measurements.
Prior
Prior
AF
Method
Method
PnP
Sparsity
RED
PnP
Sparsity
RED
Baseline
30.91
30.05
29.17
PnP
30.78
29.27
32.84
4
2.0
NPN
31.11
31.04
31.46
NPN
31.77
31.42
33.67
Bascline
27.25
26.72
27.59
PnP
24.47
23.17
24.79
8
5.0
NPN 29.88
30.02
28.84
NPN
25.72
25.32
25.34
Baseline
26.35
25.68
27.35
PnP
20.49
19.61
20.48
10.0
NPN
29.06
29.24
28.13
NPN
20.81
20.60
20.54
Table 5: MRI: PSNR (dB) for PnP vs. PnP-NPN under different priors and AFs.
Table 6: Deblurring: PSNR (dB) for PnP vs. PnP. NPN under different priors and noise levels.
In Figures 7, 8, and 9, the effect of the parameter the quality of the recovery is shown, for MRI, Deblurring, and SPC, respectively.
Fig. 10 shows the results in terms of convergence of the PnP and PnP-NPN for SPC with m/n 0.1. The color the line indicates the type of projection matrix used, orthonormal by QR (red) or designed by Eq. 3 (blue). The color shade indicates the percentage of the low-dimensional subspace p/ n, ranging from 0.1 to 0.9. The best case is with p/n. 0.1 and S designed, which was the one used in the experiments of the main paper.
Additionally, we implemented several state-of-the-an denoisers into the PnP framework. We also change the FISTA algorithm used in the paper to an alternating direction method of multipliers (ADMM) formulation. This algorithm splits the optimization problem into two subproblems: the data fidelity and the prior. We employed several -fate-of-the-ar denoisers such as Restormer [60], DnCNN [63], DnCNN-Lipschizz [43], and DRUNet [62]. We validated image deblurring application with 2.0. We used 200 iterations of the PnP-ADMM algorithm, the step size a 0.5, and the value of 7 = 0.7. Table are shown the obtained results. Here, the proposed regularization function consistently improves the baseline PnP-ADMM algorithm.
Gaussian kernel bandwidth
27

---

33.5
31.0
PnP
33.0
32.5
30.5
PSNR
FSNR
32.0
30.0
31.5
31.0
29.5
NPN
30.5
PuP
29.0
30.0
10-5
10-4
10-3
10-2
10-1
10P
10
10
Figure 7: Effect y in NPN on PSNR (dB) for MRI reconstruction, with a 1 10-4 The maximum PSNR of 33.67 dB is achieved when 7 = 6 10-3
Figure 8: Effect of ty on PSNR (dB) in deblur- ring reconstruction, with a 0.5 10-4, The maximum PSNR of 31.25 dB is achieved when 7 = 0.263.
21.00
20.75
20.50
PSNR
20.25
20.00
19.75
NPN Baseline
19.50
10
100
Figure 9: Effect of y on PSNR (dB) in SPC reconstruction, with 8 10-4, The maximum PSNR of 21.17 dB achieved when ? = 1.04.
Method
Restormer
DrCNN
DnCNN-Lipschitz
DRUNet
Sparsity Prior
Baseline
29.86
29.55
30.36
29.68
28.75
NPN
32.62
32.12
32.35
32.07
29.75
Table 7: Comparison of PnP-ADMM method for image deblurring with 2.0 using different
denoisers. The best result is highlighted in bold teal and the second- best in orange underline
A.12 Analysis on the selction of p and data adaptation
The selection of the size of S is important for the optimization (3) as the bigger value of p. the more challenging the projection estimation becomes. We show the estimation error the following table, showing that the error increases by increasing the number of rows of S.
28

---

PnP-NPN with S (QR)
PnP-NPN with S (Designed)
PnP
10
pin (Designed)
pin (QR)
0.5
0.3
10
175
200
Iteration
Figure 10: Effect of the low-dimensional subspace dimension p/n and projection matrix design on signal convergence in SPC reconstruction, with a 8 x 10
Table 8: Projection error for different values of p/n.
|Sx" GG()ll
p/r n
|Sx*]]
0.1
0.1031
0.3
0.1544
0.5
0.1821
0.7
0,2566
0.9
0.2305
Additionally, to illustrate the effect of the data-driven adaptability of S, we clarify the following. In compressed sensing, to span the null-space, we design matrix orthogonal to H via the QR decomposition, and via optimization using the regularization in Eq. (3), which promotes that the rows of S lie in the null-space of H. To analyze how the data distribution affects the design of S, we optimize four matrices, one for each dataset (MNIST, FashionMNIST,, CIFAR1O, and CelebA), cross-validating the data-driven invertibility loss Ex X S' Sx113. The results show that the QR- initialized matrix provides good invertibility across all datasets, and that optimizing S via Eq. (3) on each dataset further improves it. Thus, the data distribution of an already orthogonal matrix. Testing S on in-distribution samples yields the best invertibility (bold), while diverse training data (e.g. CIFAR10) gives the second-best (underlined), suggesting single S could generalize across datasets.
enhance the data-driven invertibility
Test/S design
QR(H) Alg
MNIST
FashionMNIST
CIFARIO
CelebA
MNIST
9.28 10
2.97 10
4.99 10
5.29 10
5.56 x
MNSHIonMNIST
1.76 10-4
7.07 x 10-3
3.11 x 10 10-5
5.05 x 10-3
5.62 x 10-5
CIFARIO
2.62 10-1
1.67 x 10
8.95
10
4.01 10
4.52 x 10
CelebA
2.39 10
1.55 x 10 9.09 10
4.68 x 10
4.15 x 10
Table 9: Cross dataset validation of data invertibility metric Ex Ex|[x S1Sx
29

---

A.13 Limitations on the selection of S
The method works when it is possible to find non-linear correlation between null space (NS) components Sx and measurements Hx, leveraging dataset with triplets (x, Hx;, Sx;). We promote this correlation by solving Eq. (3) that balances the two terms: i) ||Sx G( encourages
non-linear correlation between Sx and Hx. ii) | A TAA A| where A H [ promotes
orthogonality and ensures that S samples components from the NS of H.
There scenarios where the non-linear correlations are easier achieve, for instance, in deblurring or SR, due to the low-pass filters associated with these tasks being non-ideal, which leads to close frequency bands between the high-pass components selected by S and the Gaussian low-pass components sampled by H. In other scenarios, such as MRI or CT, the non-linear correlations between Sx and Hx more challenging to achieve due to the orthogonality of their respective rows, But, the method still works by designing $ such that leverages structural similarities of H, such adjacent frequencies in MRI neighboring angles in CT.
To illustrate the effectivvness of finding non-linear correlations between Sx and Hx, we set challenging scenario in an MRI task with H sampling only random low frequencies using ID Cartesian mask with an acceleration factor (AF) of 12. then apply an FFTSHIFT operation to this mask to produce high-frequenee frequency sampling pattern, denoted by S1, spatially distant from the frequency support of H but maintaining the same AF. We refer to this configuration as the well- separated sampling. In this case, the method fails to recover non trivial solutions, ||G (y)l 0.0108, and high relative error 1.0.
Then, we configure scenario where the non-linear correlations are easier to find, using the same H described before, we set S2 I that selects frequencies adjacent to those sampled by H. This adjacent sampling configuration corresponds precisely the setup employed in the MRI experiments presented in the main manuscript. See Figure 11. exhibits increased spectral (non- linear) correlation, resulting in small relative error of 0.325, and a norm (y)) close to the norm of the ground-truth signal [Sx|
H
S
S
Figure 11: MRI masks
SSx G
Sampling Type
|G(y)]
||Sx]
PSNR(G((), Sx) [dB]
Sx
Disjoint Sampling
1.00
0.0197
1.72
30.7
Adjacent Sampling
0.40
98.7
112
59.5
Table 10: Comparison of reconstruction metrics under different sampling strategies.
A.14 Comparison with reconstruction model
To create a fair comparison with the proposed method, we consider the network W* 1. + g" which has the same number parameters as but trained reconstruct directly the underlying signal i.e., W arg minw Ex,y ||x W(y)l. The network trained with the same number of epochs as G (300 epochs), using the Adam W optimizer with learning rate of 1e-1, We evaluate two aspects here: i) the reconstruction performance the network W and ii) incorporating into the PnPMADMM the regularization but in the form of @w (x) ||x SW [y)] We trained the
30

---

network W with the Places365 dataset. The PnP-ADMM was set with 100 iterations, using the DnCNN-Lipchitz denoiser. We set a 0.3 and y 0.2. In the following table, we report the obtained results. The results show that the proposed approach overall improves the baseline, the reconstruction of the model W, and the regularization based on this model.
Metric
Base PnP-ADMM
NPN PmP-ADMM
PnP-ADMM with Φw (y)
PSNR (dB)
22.30
23.87
21.20
19.71
SSIM
0.586
0.678
0.534
0,490
Table Comparison of PSNR (dB) and SSIM metrics
In the tables below, we show the training and execution time in seconds (s) and minutes (min), respectively. These results are obtained for a batch size of 100 images with resolution of 128 x 128 and 100 iterations the PnP-ADMM. W (reconstruction network from the measurements y) and G are trained for 300 epochs. The training times of W and G are very similar due to having the same number of parameters. Despite a modest increase of 0.20 seconds over the baseline PnP-ADMM, NPN achieves 1.57 dB PSNR gain. Moreover, compared to PnP-ADMM with éw. NPN is 0.28 seconds faster while delivering a 2.67 dB improvement in PSNR. These results show that the NPN regularization yields substantial quality gains with minimal or favorable time trade-offs.
Method
Execution time (s)
Method
Training time (min)
Base PnP-ADMM
12.48
G(y) W(y)
385
NPN PnP-ADMM
/2.
390
PnP-ADMM with dw
12.96
Table 12: Training time for the networks
Table 13: Execution times for PnP with network W and
W and G (NPN).
G (NPN).
31


Reviews:
[{"review_id": "9WDJwUi5z4", "reviewer": "Reviewer_RKuN", "summary": "The method proposes to learn a network to estimate the solution’s projection in the null space of the linear operator when solving a linear inverse problem.  It is used as an additional regularization in the optimization (along with the typical prior term) to inform the optimization process about the structure of the linear operator.", "questions": "I listed my questions above.  To reiterate, my main questions are about the construction of the low-rank S and the effectiveness of the optimization (3) to recover important subspace that is useful to describe the data distribution. I also have a question about whether in some applications the estimated G(y) would give a trivial solution.", "limitations": "yes", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths \n- The paper explains the motivations, assumptions, and the algorithm clearly. \n\n- The paper provides analyses on the convergence of utilizing the proposed regularizer in solving linear inverse problems. \n\n\nWeaknesses\n- Line 141 says the rank of S is p < (n -m), but from Sec 3.1 it seems the full basis of the null space is used.  Is S chosen to be low rank, and if so, the paper does not describe how to choose the low rank S.  In the case of DFT and random matrices, each basis vector seems equal without looking at the data distribution to be solved for.    \n\n\n- Perhaps related to the question above, how do the choice of p and the initialization of S affect the optimization (3)?  When S is initialized as a subset of the null space basis, is the optimization (3) able to get S that is composed of a different set of basis vectors? Also how does the data distribution affect the solution of S?  \n\n- I do not fully understand how the network G(y) is able to estimate Sx effectively.  My understanding is y is the measurement, which contains no information about the projection of x in the null space.  In this case, what information does G(y) use to estimate Sx?  For example, in the case of super-resolution, there will be multiple x’s that have the same measurement y but different Sx.  In this case, would G(y) simply estimate E[Sx] of the entire training set, and in the case of super-resolution if we have a large dataset will this expectation simply be 0 and we recover a typical l2 regularization in (2)?   \n\n\n- The comparison in Table 1 is mainly with other null-space methods.  How does the proposed method compare to other types of methods for solving inverse problems, e.g. [1] and [2]?  Can the proposed method be incorporated in these methods? \n\n[1] Solving Linear Inverse Problems Provably via Posterior Sampling with Latent Diffusion Models. NeurIPS 2023. \n\n[2] Solving Inverse Problems with Latent Diffusion Models via Hard Data Consistency. ICLR 2024.", "was_revised": true, "final_justification": "My concerns were (1) how S was chosen, (2) how the network G is able to estimate Sx, (3) insufficient comparison to existing methods.  The rebuttal addresses my concerns.  The authors should add the discussions in the revision.", "rebuttal": "We appreciate the reviewer’s constructive feedback, which has helped improve the paper. We also thank the reviewer for noting that the paper provides “analyses on the convergence of utilizing the proposed regularizer in solving linear inverse problems.” Below, we address each concern in detail. \n\n## Weaknesses\n\n### W1. Low rank nature of $\\mathbf{S}$\n\nWe appreciate the reviewer’s comment regarding the rank and construction of the matrix $\\mathbf{S}$. The matrix $\\mathbf{S} \\in \\mathbb{R}^{p \\times n}$ is not promoted to be low-rank; rather, it is specifically designed to have full row rank. Although we can control the number of rows, choosing $p \\leq n-m$, $\\mathbf{S}$ does not necessarily span the entire null space basis. Instead, it spans a subspace, which varies depending on the application. For example, in compressive sensing, we promote the maximum rank of $\\mathbf{A} = [\\mathbf{H}; \\mathbf{S}]$, i.e., $m+p \\leq n$.\n\nWhile we do not provide theoretical bounds for the chosen size, we evaluate its effect empirically on reconstruction quality with PnP-FISTA in Figure 8 of the supplementary material. While increasing $p$ may improve performance by capturing more null-space components, it also makes the estimation task of $\\mathrm{G}$ more challenging, as it must recover a higher-dimensional projection from measurements $\\mathbf{Hx}$. Thus, there is a tradeoff between the rank of $\\mathbf{S}$ and the network capacity required for accurate prediction.\n\nIn structured settings such as MRI with undersampled DFT matrices, we construct $\\mathbf{S}$ using the missing rows of $\\mathbf{H}$, which worked well without additional optimization. In contrast, for unstructured sensing matrices (e.g., Bernoulli in compressed sensing), optimizing $\\mathbf{S}$ via Equation (3) resulted in better convergence and reconstruction than using a fixed orthonormal basis.\n\n---\n\n### W2. Choice of $p$ and data adaptability\n\nThe size of $\\mathbf{S}$ significantly affects the optimization in Equation (3). A larger $p$ increases the difficulty of projection estimation. Using the same models as in Figure 8, we report projection errors in the table below, showing that increasing the number of rows in $\\mathbf{S}$ leads to higher estimation error:\n\n| $p/n$ | $\\Vert\\mathbf{Sx}^* - \\mathrm{G}(\\mathbf{y})\\Vert_2^2 / \\Vert\\mathbf{Sx}^*\\Vert_2^2 $ |\n|:---:|:------:|\n|0.1|**0.1031**|\n|0.5|_0.1821_|\n|0.9|0.2305|\n\nTo illustrate the effect of the data-driven adaptability of $\\mathbf{S}$, we clarify the following.  In compressed sensing, to span the null-space, we design a matrix orthogonal to $\\mathbf{H}$ via the QR decomposition, and via optimization using the regularization in Eq. (3), which promotes that the rows of $\\mathbf{S}$ lie in the null-space of $\\mathbf{H}$. To analyze how the data distribution affects the design of $\\mathbf{S}$, we optimize four matrices, one for each dataset (MNIST, FashionMNIST, CIFAR10, and CelebA), cross-validating the data-driven invertibility loss $ \\mathbb{E}_{\\mathbf{x}}\\Vert\\mathbf{x} - \\mathbf{S}^\\dagger \\mathbf{S}\\mathbf{x}\\Vert_2^2. $\nThe results show that the QR‐initialized matrix provides good invertibility across all datasets, and that optimizing $\\mathbf{S}$ via Eq. (3) on each dataset further improves it. Thus, the data distribution can enhance the data‐driven invertibility of an already orthogonal matrix. Testing $\\mathbf{S}$ on in‐distribution samples yields the $\\textbf{best}$ invertibility ($\\textbf{bold}$), while diverse training data (e.g.\\ CIFAR10) gives the _second‐best_ (_underlined_), suggesting a single $\\mathbf{S}$ could generalize across datasets.\n\n| Test \\ S design| QR | MNIST| FashionMNIST| CIFAR10| CelebA|\n|-----|:---:|:----:|:-----:|:-----:|:-----:|\n| MNIST| 9.3×10⁻⁵ | **3.0×10⁻⁵**| _5.0×10⁻⁵_| 5.3×10⁻⁵| 5.6×10⁻⁵|\n| FashionMNIST| 1.8×10⁻⁴ | 7.1×10⁻⁵ | **3.1×10⁻⁵**| _5.1×10⁻⁵_| 5.6×10⁻⁵|\n| CIFAR10| 2.6×10⁻⁴| 1.7×10⁻⁴| 9.0×10⁻⁵| **4.0×10⁻⁵**| _4.5×10⁻⁵_|\n| CelebA| 2.4×10⁻⁴ | 1.6×10⁻⁴| 9.1×10⁻⁵| _4.7×10⁻⁵_| **4.2×10⁻⁵**|\n\n---\n### W3. Training $\\mathrm{G}$ and implications\n\nWe appreciate the reviewer’s thoughtful question. Below, we clarify the role of $\\mathrm{G}$ in estimating $\\mathbf{Sx}$, particularly given that $\\mathbf{y} = \\mathbf{Hx}$ contains no explicit information about the null-space component $\\mathbf{Sx}$.\n\n1. Supervised training: The network $\\mathrm{G}$ is trained a priori, independently of the reconstruction algorithm, according to Eq (3). The dataset consists of $N$ image samples $ \\{\\mathbf{x}_i\\}$, corresponding measurements $ \\{\\mathbf{Hx}_i\\}$, and null-space projections $\\{\\mathbf{Sx}_i\\}$ for $i = 1, \\dots, N$. While $\\mathbf{Sx}$ is not directly recoverable from $\\mathbf{Hx}$, the dataset enables learning to predict $\\mathbf{Sx}$ from $\\mathbf{Hx}$ via statistical correlations. \n\n2. Regularization to avoid trivial solutions: To avoid trivial solutions $\\mathrm{G}(\\mathbf{y}) \\approx \\mathbf{0}$ or $\\mathbf{Sx = 0}, \\mathbf{x}\\neq 0$, we include orthogonality and inversion regularizers in the joint training of $\\mathbf{S}$ and $\\mathrm{G}$ (see Eq. (3)). These regularizers ensure that $\\mathbf{S}$ is meaningful and that $\\mathrm{G}$ does not collapse.\n\n3. High-frequency content recovery: In tasks like super-resolution or deblurring, high-frequency details in $\\mathbf{x}$ are lost in the measurements $\\mathbf{y} = \\mathbf{Hx}$ but preserved in $\\mathbf{Sx}$. Although this introduces ambiguity, the availability of paired samples $(\\mathbf{x}_i, \\mathbf{Hx}_i, \\mathbf{Sx}_i)$ enables $\\mathrm{G}$ to learn plausible mappings. As shown in Figure 3, the trained network $\\mathrm{G}^*$ successfully recovers high-frequency components.\n\n---\n### W4. Evaluation with diffusion models. \n\nWe thank the reviewer for this suggestion. We have integrated our proposed regularization term into two widely adopted diffusion‐model frameworks, Diffusion Posterior Sampling (DPS) [3] and DiffPIR [4]. These frameworks serve as canonical diffusion pipelines upon which newer methods build. Our regularization could likewise be incorporated into other approaches, such as latent‐space diffusion models [1][2]. Specifically, we evaluated the proposed regularization on the limited‐angle CT inverse problem, using 60 out of 180 total views (spaced every 1°) under a parallel‐beam geometry. Here $\\mathbf{H}$ corresponds to the undersampled Radon transform matrix. We construct $\\mathbf{S}$ using the missing angles of $\\mathbf{H}$, i.e., $\\mathbf{S}$ contains the remaining 120 views. We train the diffusion model for $1000$ epochs with a batch size of $4$ using the AdamW optimizer and a learning rate of 3×10⁻⁴. We used a cosine variance schedule ranging from $\\beta_1=1×10⁻⁴$ to $\\beta_T=0.02$ with $T=1000$ time steps.  The training set comprised $10,368$ slices from the LoDoPaB-CT dataset, resized to $256 \\times 256$. We evaluated the performance using 10 slices of the test set. The network $\\mathrm{G}$ corresponds to a U-Net, and it was trained for 100 epochs with a learning rate of 3×10⁻⁴  with a batch size of $4$ using the AdamW optimizer. \n\nWe incorporated the proposed NPN regularization term at line 7 of the DPS algorithm \n\n$$\\mathbf{x}^{t-1} = \\mathbf{x}^{t-1} - \\zeta_\\ell \\nabla_{\\mathbf{x}^t} \\left( \\left\\Vert \\mathbf{y} - \\mathbf{H}\\hat{\\mathbf{x}}_0 \\right\\Vert_2^2 + \\gamma \\left\\Vert\\mathrm{G}(\\mathbf{y}) - \\mathbf{S} \\hat{\\mathbf{x}}_0 \\right\\Vert_2^2 \\right)$$\nand at line 4 of the DiffPIR algorithm as follows:\n$$\\hat{\\mathbf{x}}_0^{(t)} = \\underset{\\mathbf{x}}{\\mathrm{argmin}} \\left\\Vert \\mathbf{H}\\mathbf{x} - \\mathbf{y} \\right\\Vert_2^2 + \\rho_t \\left\\Vert \\mathbf{x} - \\mathbf{x}_0^{(t)} \\right\\Vert_2^2 + \\gamma \\left\\Vert \\mathbf{S}\\mathbf{x} - \\mathrm{G}(\\mathbf{y}) \\right\\Vert_2^2$$\n\nwhere $\\gamma$ is a hyperparameter to weigh the importance of the regularization. The regularization promotes higher consistency in the sampling process. \nWe tested both DPS and DiffPIR across several $\\gamma$ values. The table below shows the obtained results. For DPS, the NPN regularization consistently improves reconstruction performance by up to 1.85 dB. For DiffPIR, it yields improvements of up to 0.61 dB (for $\\gamma = 1×10⁻⁴$). Further exploration of $\\gamma$ may yield even greater gains.\n\n| $\\gamma$ |NPN–DPS|NPN–DiffPIR|\n|-----|:----:|:--------:|\n| 0.0 (Base)| 28.22| 31.30|\n| 10⁻⁵ | 28.55  |_31.88_|\n|10⁻⁴ | 28.30  | **31.91** |\n| 10⁻³ | 28.47 | 30.53 |\n| 10⁻² | 28.78 | 29.90 |\n| 10⁻¹ |_30.06_| 28.98 |\n| 0.2| **30.07**| 28.57 |\n| 0.5| 29.90| 28.00 |\n\nAdditionally, we include the proposed regularization into the deep image prior (DIP) framework [5]. This framework uses untrained neural networks to reconstruct signals from measurements $\\mathbf{y}$. However, DIP may overfit measurement noise. Incorporating NPN regularization mitigates this by promoting consistency. The resulting optimization problem is\n\n$$\\mathrm{K}^\\ast = \\arg\\min_{\\mathrm{K}} \\Vert\\mathbf{y}-\\mathbf{H}\\mathrm{K}(\\mathbf{z})\\Vert_2^2 + \\gamma \\Vert\\mathrm{G}^\\ast(\\mathbf{y})-\\mathbf{S}\\mathrm{K}(\\mathbf{z})\\Vert_2^2, \\quad \\hat{\\mathbf{x}} = \\mathrm{K}^\\ast(\\mathbf{z}),  \\quad \\mathbf{z}\\sim\\mathcal{N}(0,1).$$\n\nThe following table shows the performance for different values of $\\gamma$, showing improvements of up to 5 dB in reconstruction quality for 10 images of the test set of the Places365 dataset. \n\n| $\\gamma$|0.0(Baseline)| 0.001| 0.01| 0.1| 1.0|\n|---|---|---|-----|----|------|\n| PSNR|17.59 ± 5.17| 18.79 ± 2.85| 21.93 ± 1.46| **22.16 ± 1.10**| 20.54 ± 2.06|\n| SSIM|0.311 ± 0.070| 0.277 ± 0.114| 0.488 ± 0.132| **0.542 ± 0.101**| 0.483 ± 0.171|\n\n[1] Solving Linear Inverse Problems Provably via Posterior Sampling with Latent Diffusion Models. NeurIPS 2023.\n\n[2] Solving Inverse Problems with Latent Diffusion Models via Hard Data Consistency. ICLR 2024.\n\n[3]  Diffusion Posterior Sampling for General Noisy Inverse Problems. ICLR 2023.\n\n[4] Denoising diffusion models for plug-and-play image restoration. CVPR 2023.\n\n[5] Deep image prior. CVPR 2018."}, {"review_id": "nOcV6PSqEh", "reviewer": "Reviewer_hHJq", "summary": "The paper proposes a new class of regularization, referred to as non-linear projections of the null space (NPN). Assuming a linear inverse problem with a sensing matrix $H$, the basic idea of NPN is to find a proper range-null space decomposition. In particular, NPN jointly searches the projection matrix and neural network restricted to the null space of $H$. The authors theoretically analyze the convergence property of the plug-and-play (PnP) algorithm with NPN and numerically demonstrate that it can converge faster than a naive PnP algorithm. Additionally, they demonstrate some image reconstruction experiments, suggesting the superiority of NPN.", "questions": "* The reviewer wonders whether the convergence of PnP-NPN is possible outside the convergence improvement zone. Although the reviewer understands that the proof of Thm. 1 requires the zone, it is preferable to analyze the convergence property of the algorithm without it, as noted above. \n\n* Although NPN successfully improves performance, it requires additional training parameters and a training process. The authors will need to compare it with existing methods in terms of computational time and/or the number of parameters.", "limitations": "Yes. However, the authors might need to add some discussions on limitations based on the additional experiments or analyses suggested by reviewers.", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 2, "significance": 3, "originality": 3, "strengths_and_weaknesses": "**Strength**\n\n* The paper proposes a new regularization technique. Although the range-null decomposition has been studied in the literature, the idea of NPN seems novel. (Quality)\n\n* The authors present both theoretical analysis and numerical demonstrations, which effectively demonstrate the superiority of NPN. (Significance)\n\n* As described in the paper, NPN is not limited to the PnP algorithm. The idea presented in this paper will appeal to a wide range of readers in machine learning and signal processing. (Significance)\n\n**Weakness**\n\n* The statement of Thm. 1 is confusing. Thm. 1 claims that the algorithm converges to a fixed point within the (finite) convergence improvement zone. However, the authors only show in the proof that the residual error decreases in the zone, which seems inadequate for showing convergence to a fixed point. A similar expression appears in l. 263: \"Theorem 1 guarantees convergence of the sequence $\\{\\mathrm x^l\\}$ to a fixed point,\"  but convergence will occur only if the zone is infinite. (Clarity, Quality)\n\n* In Sec. 5.1, the authors define the convergence improvement zone of two figures (Fig.2b and 2c) in different ways. Although both zones look to coincide, are the definitions accurate? (Clarity)\n\n* There are several typos and inaccuracies in the notations. Please check the manuscript carefully. (Clarity)\n\ne.g.) \n\n> l.159: (RNSD) -> RNSD\n\n> l.160: $\\tilde x=\\dots$ seems intuitively understandable but inaccurate.\n\n> l.172: $\\mathcal F=f_1^T\\dots$ needs $``\\dots''$.\n\n> l.194: $A=[H^T,S^T]$ should be $A=[H^T,{\\tilde S}^T]$.\n\n> l.320: \"Figure 2(b)\" should be \"2b).\"\n\n> l.876 in Appendix: \"according to ??\" should be revised.\n\n> Alg. 2-5 in Appendix: $\\frac{\\sqrt{1+4t')^2}}{2}$ should be $\\frac{\\sqrt{1+4(t')^2}}{2}$\n\n* Although PnP-NPN algorithms are given in Appendix as pseudocode, the details of the network architecture of $\\mathrm G$ are missing. It lacks the reproducibility of the paper. (Clarity)", "was_revised": true, "final_justification": "The theoretical part of the manuscript has been significantly improved thanks to the reviewers’ comments. Although the proposed method enhances convergence performance only over a finite number of iterations, the numerical results demonstrate substantial performance gains. Therefore, the reviewer considers the manuscript’s contribution to be sufficiently valuable for acceptance at the conference.", "rebuttal": "We appreciate the reviewer’s constructive feedback, which has led to a quality improvement of the paper. We also thank the reviewer for recognizing our work as \"appealing to a wide range of readers in machine learning and signal processing.\"\nBelow, we address each concern in detail.\n\n# Weaknesses\n\n## W1. Clarity on Theorem 1\n\nWe thank the reviewer for the insightful comment and agree that the fixed-point convergence claim in the original statement of Theorem 1 may be misleading. Our main contribution in this result is not to prove convergence to a unique fixed point, but rather to establish a linear residual decay rate within a bounded number of iterations, defined as the *convergence improvement zone* (CIZ):\n\n\n$$\n\\mathcal{L}_\\phi = \\{ \\ell :  \\Vert \\mathrm{N}(\\mathbf{H}\\mathbf{x}^*)\\Vert_2^2 < \\Vert\\mathbf{S}\\mathbf{x}^\\ell - \\mathbf{S}\\mathbf{x}^\\ast \\Vert_2^2 \\}\n$$\n\nSpecifically, Theorem 1 shows that the error $\\Vert\\mathbf{x}^{\\ell+1} - \\mathbf{x}^\\ast\\Vert_2^2$ decreases at a linear rate:\n\n$$\n\\rho \\triangleq (1+\\delta) \\left( \\left\\Vert \\mathbf{I} - \\alpha (\\mathbf{H}^{\\top} \\mathbf{H} + \\mathbf{S}^{\\top} \\mathbf{S}) \\right\\Vert_2^2 + (1+\\Delta_{\\mathcal{M}_D}^S) \\Vert \\mathbf{S}\\Vert_2^2\\right),\n$$\n\nwith respect to the previous iteration error $\\Vert\\mathbf{x}^{\\ell} - \\mathbf{x}^\\ast\\Vert_2^2$. Here, the dominant terms controlling the decay rate are: (i) the orthogonality between $\\mathbf{H}$ and $\\mathbf{S}$, as enforced by the optimization described in Equation (3), and (ii) the estimation error  $\\mathrm{N}(\\mathbf{H}\\mathbf{x}^{\\*})$ of the network $\\mathrm{G}^*(\\mathbf{y})$. This rate demonstrates improved stability relative to standard PnP schemes, particularly when the term $\\Vert\\mathbf{I} - \\alpha(\\mathbf{H}^{\\top}\\mathbf{H})\\Vert_2^2$ is large.\n\nTo avoid confusion, we have revised the theorem statement to explicitly remove any misleading fixed-point claim, clearly specifying that the result holds for a finite number of iterations $ L_\\phi$:\n\n**Theorem (PnP-NPN Convergence)**  \nConsider the fidelity term $g(\\tilde{\\mathbf{x}}) = \\Vert\\mathbf{y} - \\mathbf{H}\\tilde{\\mathbf{x}}\\Vert_2^2$, and assume the denoiser $D_\\sigma$ satisfies Assumption 2. Let the matrix $\\mathbf{S}$ be constructed according to Equation (3) and satisfy the RIP condition (Definition 2) with constant $\\Delta_{{M}}^S \\in [0,1)$. Suppose the network $\\mathrm{G}^*$ satisfies Assumption 3. Then, for a finite number of iterations $\\ell = 1, \\dots, L_\\phi$ within the CIZ, the residual $\\Vert\\mathbf{x}^{\\ell+1} - \\mathbf{x}^\\ast\\Vert_2^2$ decays linearly at a rate:\n\n$$\n\\rho \\triangleq (1+\\delta) \\left( \\left\\Vert \\mathbf{I} - \\alpha (\\mathbf{H}^{\\top} \\mathbf{H} + \\mathbf{S}^{\\top} \\mathbf{S}) \\right\\Vert_2^2 + (1+\\Delta_{\\mathcal{M}_D}^S) \\Vert \\mathbf{S} \\Vert_2^2 \\right).\n$$\n\n---\n## W2. Convergence improvement zones \n\nWe thank the reviewer for pointing this out. The CIZ is defined as the region in which the network estimation error satisfies $\\mathcal{L}_\\phi = \\{\\ell :\\Vert\\mathrm{N}(\\mathbf{H}\\mathbf{x}^{\\ast})\\Vert_2^2 <\\Vert\\mathbf{Sx}^\\ell - \\mathbf{Sx}^{\\ast}\\Vert_2^2\\}$. The drawn zones in Fig. 2b) and Fig. 2c) should match as the algorithm, as Fig. 2c) is the theoretical definition of the CIZ, and Fig. 2b) is the empirical validation of the improvement in the convergence ratio. \n\nIn Fig. 2b, we illustrate the behavior of the algorithm using the empirical convergence ratio $\\Vert\\mathbf{x}^{\\ell+1} - \\mathbf{x}^\\ast\\Vert_2^2 / \\Vert\\mathbf{x}^\\ell - \\mathbf{x}^\\ast\\Vert_2^2$. We acknowledge an error in the graphical annotation of the convergence zone: for instance, in the blue curve, the proposed regularized algorithm with $\\mathbf{S}$ designed in Eq. (3), the zone should terminate around iteration $\\ell = 60$, where the black curve (baseline algorithm) drops below the blue curve. Similarly, for the red curve, $\\mathbf{S}$  (QR), the zone should go until around $\\ell = 120$.\n\nIn Fig. 2c, the CIZ is defined based on the inequality $\\mathcal{L}_\\phi = \\{\\ell :\\Vert\\mathrm{N}(\\mathbf{H}\\mathbf{x}^\\ast)\\Vert_2^2 <\\Vert\\mathbf{Sx}^\\ell - \\mathbf{Sx}^{\\ast}\\Vert_2^2\\}$, indicated by the intersection of the solid and dotted lines. The key observation across both figures is that the convergence of the algorithm is enhanced in the regime where the estimation error of the network is smaller than the projection error, indicating that the learned regularization is providing useful information to guide the recovery process. With the correction mentioned earlier, the zones in both figures align more closely, validating the developed theory.\n\n---\n## W3. Typos in the manuscript\n\nWe appreciate the reviewer's careful reading of the paper. We will fix the raised typos in the final version of the paper. \n\n---\n## W4. Network architecture details\n\nWe appreciate the reviewer’s concern regarding clarity and reproducibility. We would like to kindly clarify that Section A.5 of the Supplementary Material contains the architecture details of the network $\\mathrm{G}$ used in all experiments. To improve clarity, we will explicitly state in the revised Supplementary that this section refers to $\\mathrm{G}$. However, these are the model specifications used in the experiments.\n\n **SPC**: For the single‐pixel camera experiments, we employ a ConvNeXt‐inspired backbone [1]. The network has five ConvNeXt blocks, each comprising two successive $3\\times3$ convolutions with ReLU activations. A final $3\\times3$ convolutional layer projects the 128 features back to one channel. The output feature map is flattened and passed through a linear module to match the measurement dimensionality. We use cosine‐based positional encoding with a dimension of 256, and set the number of blocks to 5 and the base channel width to 128. \n\n**MRI, Deblurring and CT**: To train $\\mathrm{G}$ for the MRI and Deblurring experiments, we use a U-Net architecture [2] with three downscaling and three upscaling modules. Each module consists of two consecutive Conv → ReLU blocks. The downscaling path uses max pooling for spatial reduction, starting with 128 filters and increasing up to 1,024 at the bottleneck. The upscaling path performs nearest-neighbor interpolation before each module, progressively reducing the number of filters. A final 2D convolutional layer without activation produces the output. Skip connections link corresponding layers in the encoder and decoder.\n\n[1] Zhuang Liu, et al. A ConvNet for the 2020s. In CVPR 2022.\n\n[2] Olaf Ronneberger, et al. U-net: Convolutional networks for biomedical image segmentation. In MICCAI 2015.\n\n---\n# Questions \n\n## Q1. Convergence without an improvement zone\n\nWe thank the reviewer for raising this important point. The convergence of PnP-NPN outside the CIZ can indeed be studied under the general theoretical frameworks developed in the PnP literature. For instance, [1, Theorem 1] for PnP-projected gradient descent and [2, Theorem 1] for PnP-ADMM provide convergence guarantees under assumptions that are also satisfied in our setting, such as boundedness of the denoiser and Lipschitz continuity of the data fidelity gradient.\n\nHowever, our inclusion of the CIZ serves to highlight a specific regime in which the proposed NPN regularization actively contributes meaningful signal information — particularly in the form of structured null-space projections. Within this zone, we can rigorously characterize an improved convergence rate due to the alignment between the learned prior and the null-space component of the signal.\n\nWhile convergence outside this zone is still covered by general PnP theory, our main contribution lies in identifying and exploiting this beneficial regime where convergence is not only guaranteed but accelerated due to the informative nature of the learned null-space prior.\n\n[1] Provable preconditioned plug-and-play approach for compressed sensing MRI reconstruction. IEEE TCI, 2024. \n\n[2] Plug-and-play ADMM for image restoration: Fixed-point convergence and applications. IEEE TCI, 2016.\n\n---\n\n## Q2. Comparison in terms of computational time and number of parameters\n\nThis is a very nice question. To create a fair comparison with the proposed method, we consider the network $ W^* : \\mathbb{R}^m \\to \\mathbb{R}^n$  which has the same number of parameters as $G$ but is trained to reconstruct directly the underlying signal, i.e. $ W^* = \\arg\\min_{W}\\ \\mathbb{E}_{\\mathbf{x},\\mathbf{y}}\\big\\Vert\\mathbf{x} - W(\\mathbf{y})\\big\\Vert_2^2.$  \n\nThe network was trained for 300 epochs (same as $G$) using the AdamW optimizer with learning rate $1\\times10^{-4}$. We evaluate two aspects:\n\n1. The reconstruction performance of the network $W$.  \n2. Incorporation into PnP-ADMM of the regularization in the form  $$\n   \\phi_{W}(\\mathbf{x}) = \\big\\Vert\\mathbf{S}\\mathbf{x} - \\mathbf{S}W^*(\\mathbf{y})\\big\\Vert_2^2.\n   $$\n\nWe trained $W$ on the Places365 dataset. PnP-ADMM was run for 100 iterations with the DnCNN-Lipschitz denoiser, setting $\\alpha = 0.3$ and $\\gamma = 0.2$.  The results show that the proposed approach overall improves the baseline,  the reconstruction of the model $\\mathrm{W}$, and the regularization based on this model.\n\n| Metric| Base PnP-ADMM | NPN PnP-ADMM | PnP-ADMM with $\\phi_{W}$ | $W^*(\\mathbf{y})$ |\n|--------------:|--------------:|-------------:|-----------------------:|------------------:|\n|PSNR|_22.30_ |**23.87** |21.20 |19.71 |\n|SSIM|_0.586_|**0.678** |0.534|0.490|"}, {"review_id": "MwH0HXNNRI", "reviewer": "Reviewer_PAt7", "summary": "The paper introduces a method for solving imaging inverse problems by leveraging the structure of the data in the null space of the measurement operator. The authors propose a novel regularization approach that uses a neural network trained to estimate the null space component from measurements in a supervised manner. Reconstruction is then performed via a Plug-and-Play (PnP) scheme that uses the gradient of the learned regularization.", "questions": "### Suggestions\n\n**Theory**\n\n- Clearly state whether the convergence result of Theorem 1 applies to the infinite sequence or only a finite number of iterations.\n\n- Provide practical guidance on how to choose step size $\\alpha$ and offer examples or justification for when the operator $S$ satisfies the required conditions (e.g., low spectral norm, low RIC $\\Delta^S_{M_D}$).\n\n- Include theoretical or empirical analysis to demonstrate that the convergence improvement zone is non-empty and reachable in practice.\n\n**Experiments**\n\n- Improve baseline comparisons and include more datasets and problems to give more information to the reader.\n\n- Provide results showing the effect of different PnP denoisers, step sizes, and regularization parameters ($\\lambda_1$, $\\lambda_2$) and discuss how sensitive the method is to these hyperparameters in practice.", "limitations": "The authors discussed limitations.", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 4, "significance": 2, "originality": 3, "strengths_and_weaknesses": "### Strengths\n\n- The proposed method leverages the structure of the null space, a promising and active direction in inverse problems.\n\n- The paper is clearly written and well-organized, making the technical content accessible.\n\n### Weaknesses\n\n- The theoretical analysis and its practical relevance are somewhat unclear:\n    * In Theorem 1, the authors claim that the algorithm converges to a fixed point, but the convergence analysis is limited to a possibly finite number of iterations.\n    * Lines 239–251 offer informal justification for the assumptions required for convergence. However, the paper lacks both practical and theoretical evidence that these assumptions hold in realistic settings. For example, what would be an appropriate step size $\\alpha$? What kind of operator satisfies the condition that the spectral norm of $S$ and its RIC $\\Delta^S_{M_D}$ are sufficiently small?\n    * Since the \"convergence\" behavior only holds within the \"convergence improvement zone,\" it would be helpful to know whether this zone is non-empty and practically reachable. Can the authors provide empirical or theoretical evidence that iterates typically enter this region?\n\n- The experimental section could be significantly strengthened:\n    * The chosen baseline appears weak. More recent and competitive methods are readily available (see https://deepinv.github.io/deepinv/ for a list of algorithms and networks).\n    * It would benefit from additional datasets and a wider variety of inverse problems and/or parameterizations.\n    * The sensitivity of the method to key hyperparameters—such as step size, PnP denoiser choice, and regularization weights $\\lambda_1$ and $\\lambda_2$—is not thoroughly explored. From a practical standpoint, it would be valuable to understand how robust the method is to these choices.", "was_revised": true, "final_justification": "I upgraded my rating after the rebuttal and discussion with the authors, which helped clarify the practical interest of the theory. I recommend acceptance of the paper, though the authors should include all the details given in the discussion to the paper to make it clearer.", "rebuttal": "We appreciate the reviewer’s constructive feedback, which has helped improve the paper. We also thank the reviewer for noting that the proposed method is “a promising and active direction in inverse problems”. Below, we address each concern in detail.  \n\n# Weaknesses\n---\n## W1. Clarity on Theorem 1\nWe thank the reviewer for the insightful comment and agree that the fixed-point convergence claim in the original statement of Theorem 1 may be misleading. Our main contribution in this result is not to prove convergence to a unique fixed point, but rather to establish a linear residual decay rate within a bounded number of iterations, defined as the *convergence improvement zone* (CIZ)  Specifically, Theorem 1 shows that the error $\\Vert\\mathbf{x}^{\\ell+1} - \\mathbf{x}^\\ast\\Vert_2^2$ \ndecreases at a linear rate:\n$$\\rho \\triangleq (1+\\delta) \\left( \\left\\Vert \\mathbf{I} - \\alpha (\\mathbf{H}^{\\top} \\mathbf{H} + \\mathbf{S}^{\\top} \\mathbf{S}) \\right\\Vert_2^2 + (1+\\Delta_{\\mathcal{M}_D}^S) \\Vert \\mathbf{S}\\Vert_2^2\\right),$$\nwith respect to the previous iteration error $\\Vert\\mathbf{x}^{\\ell} - \\mathbf{x}^\\ast\\Vert$. Here, the dominant terms controlling the decay rate are: (i) the orthogonality between $\\mathbf{H}$ and $\\mathbf{S}$, as enforced by the optimization described in Equation (3), and (ii) the estimation error  $\\mathrm{N}(\\mathbf{H}\\mathbf{x}^{\\*})$ of the network $\\mathrm{G}^*(\\mathbf{y})$. \n\n**Theorem (PnP-NPN Convergence)**  \nConsider the fidelity term $g(\\tilde{\\mathbf{x}}) = \\Vert\\mathbf{y} - \\mathbf{H}\\tilde{\\mathbf{x}}\\Vert_2^2$, and assume the denoiser $D_\\sigma$ satisfies Assumption 2. Let the matrix $\\mathbf{S}$ be constructed according to Equation (3) and satisfy the RIP condition (Definition 2) with constant $\\Delta_{{M}}^S \\in [0,1)$. Suppose the network $\\mathrm{G}^*$ satisfies Assumption 3. Then, for a finite number of iterations $\\ell = 1, \\dots, L_\\phi$ within the CIZ, the residual $\\Vert\\mathbf{x}^{\\ell+1} - \\mathbf{x}^\\ast\\Vert_2^2$ decays linearly at a rate:\n$$\n\\rho \\triangleq (1+\\delta) \\left( \\left\\Vert \\mathbf{I} - \\alpha (\\mathbf{H}^{\\top} \\mathbf{H} + \\mathbf{S}^{\\top} \\mathbf{S}) \\right\\Vert_2^2 + (1+\\Delta_{\\mathcal{M}_D}^S) \\Vert \\mathbf{S} \\Vert_2^2 \\right).$$\n\n---\n## W2. Assumptions in theoretical derivations\nThanks for these interesting questions. The step size $\\alpha$ was chosen to $\\alpha = 1/\\Vert\\nabla^2 g(\\mathbf{x})\\Vert_2^2$ following theory on gradient descent [1]. Additionally, in the compressed sensing theory, a wide range of studies have been performed on computing RIC for random Bernoulli or Gaussian matrices (the ones used in CS experiments) [2], proving that RIC is around 0.1-0.3 with high probability. \n\n[1] Notes for Optimization Algorithms Spring 2023 (2023).\n\n[2] Compressive sensing and structured random matrices (2010).\n\n---\n## W3. Convergence Improvement Zone (CIZ)\nEmpirical evidence provided in Fig. 2 confirms that the CIZ is non-empty. This zone is important because, within it, the learned regularizer contributes meaningful reconstruction information, enhancing convergence and accuracy. For clarity, we present numerical values extracted from Fig. 2(c) in the Table below. The network estimation error was $\\Vert\\mathrm{N}(\\mathbf{H}\\mathbf{x}^\\ast)\\Vert_2^2=  0.0584$. The CIZ, therefore, corresponds to iterations where \\$\\Vert\\mathbf{Sx}^\\ell - \\mathbf{Sx}^\\ast\\Vert_2^2 > 0.0584\\$. Representative iterations from the algorithm were selected, evaluating projection error $||\\mathbf{Sx}^\\ast - \\mathbf{Sx}^\\ell ||_2^2$ for NPN and baseline methods. As shown, when projection error exceeds network estimation error, NPN demonstrates consistently lower convergence ratios and higher PSNR values compared to baseline, indicating superior convergence performance.\n\n|\\$\\ell\\$ | NPN $\\Vert\\mathbf{Sx}^\\* - \\mathbf{Sx}^\\ell\\Vert_2^2\\$ | Base $\\Vert\\mathbf{Sx}^\\* - \\mathbf{Sx}^\\ell\\Vert_2^2\\$ | NPN PSNR | Base PSNR  | Base Convergence Ratio|NPN Convergence Ratio|\n| -------- | ------- | -----|---|---|---|------|\n|5|5.124|5.35|5.96|5.82|0.9970|0.9938|\n|50|0.814|2.80|12.53|8.13|0.9914|0.9717|\n|75|0.077|1.48|18.65|10.32|0.9885|0.9822|\n|150|0.040|0.09|20.82|18.08|0.9927|0.9995|\n\nThe following table presents results for the deblurring task with a Gaussian kernel ($\\sigma = 2.0$) using PnP-ADMM and a DnCNN denoiser. In this case, since the estimation error of the network is quite low, $\\Vert\\mathrm{N}(\\mathbf{H}\\mathbf{x}^\\ast)\\Vert_2^2 = 1.37×10⁻⁴$,  the CIZ covers all the algorithm iterations. Thus, the NPN method consistently achieves better convergence ratios and yields improved reconstructions in all iterations. \n\n $\\ell$ | NPN $\\Vert\\mathbf{Sx}^\\ast - \\mathbf{Sx}^\\ell\\Vert_2^2$| Base $\\Vert\\mathbf{Sx}^\\* - \\mathbf{Sx}^\\ell\\Vert_2^2$| NPN PSNR | Base PSNR  | Base Convergence Ratio| NPN Convergence Ratio |\n|-----|----|-----|----|-----|---|----|\n| 5|5.51×10⁻⁴|8.44×10⁻⁴| 30.84|29.56|0.9932|0.9774|\n| 25|3.81×10⁻⁴|8.70×10⁻⁴| 32.00|29.54|1.0007|0.9983|\n| 50|3.57×10⁻⁴|8.79×10⁻⁴| 32.21| 29.47|1.0000|0.9995|\n| 75|3.49×10⁻⁴|8.66×10⁻⁴| 32.32| 29.49|0.9999|0.9995|\n\nBased on these results, we conclude that the CIZ is not empty and is achievable in practice. Additionally, analysis demonstrates that an effective estimation of the network $\\mathrm{G}$ provides larger improvement margins.\n\n---\n## W4. Additional reconstruction methods\nWe appreciate your suggestion. We consider expanding the application of the proposed method to two new reconstruction frameworks.\n1. Deep image prior (DIP) [1] reconstruction framework for the deblurring task using a new dataset, Places365.\n2. Two diffusion models: reconstruction frameworks, DPS [2] and DiffPIR [3], for limited-angle CT.\n\n### Deep Image Prior\nWe study the deblurring task using the Places365 dataset, with 28,000 training and 7,000 testing images, all resized to 128×128. The DIP framework uses untrained neural networks to reconstruct signals from measurements $\\mathbf{y}$. However, DIP may overfit measurement noise. Incorporating NPN regularization mitigates this by promoting consistency. The resulting optimization problem is:\n$$\\mathrm{K}^\\ast = \\arg\\min_{\\mathrm{K}} \\Vert\\mathbf{y}-\\mathbf{H}\\mathrm{K}(\\mathbf{z})\\Vert_2^2 + \\gamma \\Vert\\mathrm{G}^\\ast(\\mathbf{y})-\\mathbf{S}\\mathrm{K}(\\mathbf{z})\\Vert_2^2, \\quad \\hat{\\mathbf{x}} = \\mathrm{K}^\\ast(\\mathbf{z}),  \\quad \\mathbf{z}\\sim\\mathcal{N}(0,1).$$\nThe following table shows the performance for different values of $\\gamma$, showing improvements of up to 5 dB in reconstruction quality for 10 images of the test set of the Places365 dataset. \n|$\\gamma$|0.0(Base)|0.001|0.01|0.1|0.5|1.0|\n|----|----|--|----|---|---|---|\n|PSNR|17.59±5.2|18.79±2.9|_21.93±1.5_| **22.16±1.1**| 21.19±1.9|20.54±2.1|\n|SSIM|0.31±0.07| 0.28±0.11| 0.49±0.13| **0.54±0.10**|_0.5±0.14_|0.48±0.17|\n\n---\n### Diffusion Models (DMs)\n\nWe consider two types of DMs: DPS [2] and DiffPIR [3]. These frameworks have become standard in solving imaging inverse problems with diffusion priors. In DPS, a measurement consistency function is used to condition the generative process. We incorporate the NPN regularization into the posterior sampling step. For every time-step $t$, we modify the posterior sampling step as:\n$$ \\mathbf{x}^{t-1} ={\\mathbf{x}^{t-1}} - \\zeta_\\ell \\nabla_{\\mathbf{x}^t} \\left(\\left\\Vert \\mathbf{y} - \\mathbf{H}\\hat{\\mathbf{x}}_0 \\right\\Vert_2^2 + \\gamma \\Vert\\mathrm{G}(\\mathbf{y}) -\\mathbf{S}\\hat{\\mathbf{x}}_0\\Vert_2^2 \\right)$$\nIn DiffPIR, we modify the data proximal subproblem with:\n$$\\hat{\\mathbf{x}}_0^{(t)} =\\underset{\\mathbf{x}}{\\mathrm{argmin}} \\left\\Vert \\mathbf{H}\\mathbf{x} - \\mathbf{y} \\right\\Vert_2^2 + \\rho_t \\left\\Vert \\mathbf{x} - \\mathbf{x}_0^{(t)} \\right\\Vert_2^2 + \\gamma \\left\\Vert \\mathbf{S}\\mathbf{x} - \\mathrm{G}(\\mathbf{y}) \\right\\Vert_2^2$$\n\nWe train the DM for 1000 epochs, batch size 4, using AdamW with learning rate $3×10⁻⁴$. Noise schedule: cosine from $\\beta_1=10⁻⁴$ to $\\beta_T=0.02,  T=1000 $. Training data: 10,368 slices from the LoDoPaB-CT dataset, resized to 256x256. Evaluation on 10 test slices.\n| $ \\gamma$| NPN–DPS | NPN–DiffPIR |\n|-----|-----|----|\n|0.0(Baseline)|28.22|_31.30_|\n|10⁻⁴|28.30|**31.91**|\n|10⁻³|28.47|30.53|\n|0.1|_30.06_|28.98|\n|0.2|**30.07**|28.57|\n|0.5|29.90|28.0|\n\n[1] Deep image prior. CVPR 2018.\n\n[2] Diffusion Posterior Sampling for General Noisy Inverse Problems. ICLR 2023.\n\n[3] Denoising diffusion models for plug-and-play image restoration. CVPR 2023.\n\n---\n## W5. Additional PnP denoisers\n\nWe thank the reviewer for the comments. Section A.8 of the supplementary material presents results for various $\\gamma$ values across SPC, MRI, and deblurring tasks, showing consistent improvements over the baseline and demonstrating robustness to the regularization parameter. We also replaced FISTA with ADMM in the PnP framework, separating data fidelity and prior subproblems, and evaluated several denoisers: Restormer [1], DnCNN [2], DnCNN-Lipschitz [3], and DRUNet [4]. For deblurring with a Gaussian kernel ($\\sigma = 2.0$), we ran 200 ADMM iterations with step size $\\alpha = 0.5$ and $\\gamma = 0.7$. Results show that our regularization improves performance across all tested denoisers.\n| Method| Restormer | DnCNN | DnCNN-Lipschitz | DRUNet | Sparsity|\n|---|---|----|-----|---|--|\n| Baseline| 29.86| 29.55 | 30.36| 29.68  | 28.75|\n| NPN| **32.62** | 32.12 | _32.35_| 32.07| 29.75|\n\nHyperparameters λ₁ = 0.001 and λ₂ = 0.01, which were found to be optimal, did not significantly affect performance. The most influential hyperparameter is $\\gamma$, the regularizer weight, as demonstrated in DMs, PnP, and DIP reconstruction frameworks.\n\n[1] Restormer: Efficient Transformer for High-Resolution Image Restoration. CVPR 2022.\n\n[2] Beyond a Gaussian denoiser: Residual learning of deep CNN for image denoising. IEEE TIP 2017.\n\n[3] Plug-and-Play Methods Provably Converge with Properly Trained Denoisers. ICML 2019.\n\n[4] Plug-and-Play Image Restoration with Deep Denoiser Prior. IEEE TPAMI 2021."}]

Author rebuttal:

Meta-review: All three reviewers who provided a full review recommend acceptance (two borderline). The paper proposes a novel regularization technique for imaging inverse problems by focusing on the null space of the sensing operator. The reviewers praised the clear motivation, theoretical analysis, and the novelty of the approach. The main weaknesses identified were related to the clarity of some theoretical aspects and the need for more extensive comparisons, but the authors' rebuttals successfully addressed these concerns.
Decision: Accept (poster)