Society for Industrial and Applied Mathematics: SIAM Journal on Imaging Sciences: Table of Contents

On the Convergence of the Iterative Regularization Method Assisted by the Graph Laplacian with Early Stopping

Harshit Bajpai — 2026-04-02T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 643-676, June 2026.
Abstract.We present a data-assisted iterative regularization method for solving ill-posed inverse problems. The proposed approach, termed IRMGL+[math], integrates classical iterative techniques with a data-driven regularization term realized through an iteratively updated graph Laplacian. Our method commences by computing a preliminary solution using any suitable reconstruction method, which then serves as the basis for constructing the initial graph Laplacian. The solution is subsequently refined through an iterative process, where the graph Laplacian is simultaneously recalibrated at each step to effectively capture the evolving structure of the solution. A key innovation of this work lies in the formulation of this iterative scheme and the rigorous justification of the classical discrepancy principle as a reliable early stopping criterion specifically tailored to the proposed method. Under standard assumptions, we establish stability and convergence results for the scheme when the discrepancy principle is applied. Furthermore, we demonstrate the robustness and effectiveness of our method through numerical experiments utilizing four distinct initial reconstructors [math]: the adjoint operator (Adj), filtered back projection, total variation denoising, and standard Tikhonov regularization. It is observed that IRMGL + Adj demonstrates a distinct advantage over the other initializers, producing a robust and stable approximate solution directly from a basic initial reconstruction.

Normalized Radon Cumulative Distribution Transforms for Invariance and Robustness in Optimal Transport Based Image Classification

Matthias Beckmann — 2026-04-06T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 677-709, June 2026.
Abstract.The Radon cumulative distribution transform (R-CDT) is an easy-to-compute feature extractor that facilitates image classification tasks especially in the small data regime. It is closely related to the sliced Wasserstein distance and provably guarantees the linear separability of image classes that emerge from translations or scalings. In many real-world applications, like the recognition of watermarks in filigranology, however, the data is subject to general affine transformations originating from the measurement process. To overcome this issue, we recently introduced the so-called max-normalized R-CDT that only requires elementary operations and guarantees the separability under arbitrary affine transformations. The aim of this paper is to continue our study of the max-normalized R-CDT especially with respect to its robustness against nonaffine image deformations. Our sensitivity analysis shows that its separability properties are stable provided the Wasserstein-infinity distance between the samples can be controlled. Since the Wasserstein-infinity distance only allows small local image deformations, we moreover introduce a mean-normalized version of the R-CDT. In this case, robustness relates to the Wasserstein-2 distance and also covers image deformations caused by impulsive noise, for instance. Our theoretical results are supported by numerical experiments showing the effectiveness of our novel feature extractors as well as their robustness against local nonaffine deformations and impulsive noise.

Shape-Informed Graph Neural Networks and Data Assimilation: Application to Velocity and Pressure Reconstruction in Aortic Blood Flow

Francesco Romor — 2026-04-22T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 710-751, June 2026.
Abstract.Computational hemodynamics can enhance image-based diagnosis and provide complementary insights to predict, understand, and monitor treatments. The high computational costs and the complexity associated with handling patient-specific settings remain a major challenge toward clinical applications. In this work, we propose a novel robust shape registration method for nonparametric aortic geometries, describing different applications for projection-based reduced-order modeling for the training of graph neural networks and for data assimilation. The registration approach is based on ResNet-LDDMM, trained with a dataset of synthetic shapes, generated from real ones with statistical shape modeling. The optimization is tailored to surface meshes and does not rely on a priori assumptions on domain parameterization. We employ a multigrid strategy during the training phase that allows handling realistic mesh sizes. The registration enables the definition of geometric encoding of different blood flow solutions on a single reference shape, as well as the design of projection-based reduced-order models. We use this geometrical encoding to improve the training of graph neural networks and to present potential applications in data assimilation problems, combined with a generalized parameterized-Background Data-Weak formulation. As a particular example of data assimilation problem, we address the reconstruction of velocity fields and wall shear stresses, as well as the estimation of pressure fields and pressure-related biomarkers, such as the pressure drop, from low-resolution velocity observations. We show various numerical tests based on synthetic data, comparing the proposed strategies with state-of-the-art estimators.

Rethinking Coupled Tensor Analysis for Hyperspectral Superresolution: Recoverable Modeling under Endmember Variability

Meng Ding — 2026-04-23T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 752-781, June 2026.
Abstract.This work revisits the hyperspectral superresolution (HSR) problem, i.e., fusing a pair of spatially co-registered hyperspectral (HSI) and multispectral (MSI) images to recover a superresolution image (SRI) that enhances the spatial resolution of the HSI. Coupled tensor decomposition–based methods have gained traction in this domain, offering recoverability guarantees under various assumptions. Existing models such as canonical polyadic decomposition (CPD) and Tucker decomposition provide strong expressive power but lack physical interpretability. The block-term decomposition model with rank-[math] terms (the LL1 model) yields interpretable factors under the linear mixture model (LMM) of spectral images, but LMM assumptions are often violated in practice—primarily due to nonlinear effects such as endmember variability (EV). To address this issue, we propose representing spectral images using a more flexible block-term tensor model with rank-[math] terms (the LMN model). This modeling choice retains interpretability, subsumes CPD, Tucker, and LL1 as special cases, and robustly accounts for non-ideal effects such as EV, offering a balanced trade-off between expressiveness and interpretability for HSR. Importantly, under the LMN model for HSI and MSI, recoverability of the SRI can still be established under proper conditions—providing strong theoretical support. Extensive experiments on synthetic and real datasets demonstrate the effectiveness and robustness of the proposed method.

A Primal-Dual Algorithm for Image Reconstruction with Input-Convex Neural Network Regularizers

Matthias J. Ehrhardt — 2026-04-23T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 782-806, June 2026.
Abstract.We address the optimization problem in a data-driven variational reconstruction framework, where the regularizer is parametrized by an input-convex neural network. While gradient-based methods are commonly used to solve such problems, they struggle to effectively handle nonsmooth problems, which often leads to slow convergence. Moreover, the nested structure of the neural network complicates the application of standard nonsmooth optimization techniques, such as proximal algorithms. To overcome these challenges, we reformulate the problem and eliminate the network’s nested structure. By relating this reformulation to epigraphical projections of the activation functions, we transform the problem into a convex optimization problem that can be efficiently solved using a primal-dual algorithm. We also prove that this reformulation is equivalent to the original variational problem. Through experiments on several imaging tasks, we show that the proposed approach not only outperforms subgradient methods and even accelerated methods in the smooth setting but also facilitates the training of the regularizer itself.

Spectral-Spatial Extraction through Layered Tensor Decomposition for Hyperspectral Anomaly Detection

Quan Yu — 2026-04-27T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 807-838, June 2026.
Abstract.Low rank tensor representation (LRTR) methods are very useful for hyperspectral anomaly detection (HAD). However, existing LRTR methods often overlook spectral anomaly and rely on computationally expensive large-scale matrix singular value decomposition. To overcome these limitations, we propose the following highly efficient layered tensor decomposition (LTD) framework that simultaneously optimizes two key components within a unified model: layer 1, which reduces spectral redundancy and extracts spectral anomaly, and layer 2, which captures spatial low rank features and extracts spatial anomaly. The resulting spectral and spatial anomaly maps are then integrated to achieve a robust final detection result. An iterative algorithm based on proximal alternating minimization is developed to solve the proposed LTD model, with convergence guarantees provided. Moreover, we introduce a rank reduction strategy with validation mechanism that adaptively reduces data size while preventing excessive reduction. Theoretically, we rigorously establish the equivalence between the tensor tubal rank and tensor group sparsity regularization (TGSR) and, under mild conditions, demonstrate that the relaxed formulation of TGSR shares the same global minimizers and optimal values as its original counterpart. Experimental results on the Airport-Beach-Urban and MVTec datasets demonstrate that our approach outperforms state-of-the-art methods.

Smooth Optimization Using Global and Local Low-Rank Regularizers

Rodrigo A. Lobos — 2026-04-30T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 839-878, June 2026.
Abstract.Many inverse problems and signal processing problems involve low-rank regularizers based on the nuclear norm. Commonly, proximal gradient methods (PGMs) are adopted to solve this type of nonsmooth problems, as they can offer fast and guaranteed convergence. However, PGM methods cannot be simply applied in settings where low-rank models are imposed locally on overlapping patches; therefore, heuristic approaches have been proposed that lack convergence guarantees. In this work we propose to replace the nuclear norm with a smooth approximation in which a Huber-type function is applied to each singular value. By providing a theoretical framework based on singular value function theory, we show that important properties can be established for the proposed regularizer, such as convexity, differentiability, and Lipschitz continuity of the gradient. Moreover, we provide a closed-form expression for the regularizer gradient, enabling the use of standard iterative gradient-based optimization algorithms (e.g., the nonlinear conjugate gradient) that can easily address the case of overlapping patches and have well-known convergence guarantees. In addition, we provide a novel step-size selection strategy based on a quadratic majorizer of the line search function that leverages the Huber characteristics of the proposed regularizer. Finally, we assess the proposed optimization framework by providing empirical results in dynamic magnetic resonance imaging (MRI) reconstruction in the context of local low-rank models with overlapping patches.

Total Normal Curvature Regularization and Its Minimization for Surface and Image Smoothing

Tianle Lu — 2026-04-30T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 879-912, June 2026.
Abstract.We introduce a novel formulation for curvature regularization by penalizing normal curvatures from multiple directions. This total normal curvature regularization is capable of producing solutions with sharp edges and precise isotropic properties. To tackle the resulting high-order nonlinear optimization problem, we reformulate it as the task of finding the steady-state solution of a time-dependent partial differential equation (PDE) system. Time discretization is achieved through operator splitting, where each subproblem at the fractional steps either has a closed-form solution or can be efficiently solved using advanced algorithms. Our method circumvents the need for complex parameter tuning and demonstrates robustness to parameter choices. The efficiency and effectiveness of our approach have been rigorously validated in the context of surface and image smoothing problems.

Multilevel Bregman Proximal Gradient Descent

Yara Elshiaty — 2026-05-05T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 913-942, June 2026.
Abstract.We present the Multilevel Bregman Proximal Gradient Descent (ML-BPGD) method, a novel multilevel optimization framework tailored to constrained convex problems with relative Lipschitz smoothness. Our approach extends the classical multilevel optimization framework (MGOPT) to Bregman-based geometries and constrained domains. We provide a rigorous analysis of ML-BPGD for multiple levels and establish well-posedness together with a global linear convergence rate. We demonstrate the effectiveness of ML-BPGD in the context of image reconstruction and validate its performance through numerical experiments.

Ground Penetrating Synthetic Aperture Imaging of Subsurface Targets

Pedro González-Rodríguez — 2026-05-05T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 943-974, June 2026.
Abstract.Motivated by the use of unmanned aerial vehicles (UAVs) equipped with ground penetrating synthetic aperture radar (GP-SAR) for imaging buried landmines, we study imaging of subsurface targets. A key challenge in this problem is the removal of so-called ground bounce signals which are reflections from the air-soil interface. When ground bounce signals can be removed exactly from GP-SAR measurements, one can readily image subsurface targets using traditional synthetic aperture imaging methods such as Kirchhoff migration (KM). However, ground bounce signals dominate measurements, carry no useful information about subsurface targets, and disrupt their imaging. In this paper we present an asymptotic theory for this problem assuming that the elevation of the UAV flight path is the largest length scale in the problem. Using this asymptotic theory, we introduce data-driven methods for (i) estimating the ground permittivity, and (ii) approximately removing ground bounce signals from GP-SAR. We first establish these methods for the simple case of a flat air-soil interface and then apply them to an unknown rough air-soil interface. We validate this asymptotic theory using numerical simulations.

Robust One-Bit Compressed Sensing via Continuous Sign Approximation and Hybrid Ordinary-Welsch Function

Zihao He — 2026-05-08T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 975-1014, June 2026.
Abstract.One-bit quantized measurements are commonly encountered in resource-constrained imaging applications such as medical imaging and satellite remote sensing. Although stable signal reconstruction from noise-free one-bit measurements has been widely demonstrated, noise in practical data acquisition systems often leads to unknown sign flips, significantly increasing the difficulty of accurate recovery. While numerous robust methods have been developed for signal reconstruction from noisy one-bit measurements, most of them necessitate prior knowledge of the number of sign flips, which is typically unavailable in practice. In this paper, we propose two novel optimization models based on continuous approximation functions (CAFs) of the sign function. The first model employs the quadratic loss function, specifically designed for scenarios with low sign flipping ratios, while the second incorporates a Hybrid Ordinary-Welsch (HOW) loss function, yielding enhanced robustness under higher flipping ratios without requiring prior knowledge. Then, we propose two algorithmic frameworks, OBCSA_CAF and OBCSA_CAF_HOW, to solve the respective models, and prove that both algorithms generate convergent sequences of objective function values and subsequences that, with high probability, converge to specific accumulation points. Extensive simulations demonstrate the superior performance of our methods over the state-of-the-art algorithms in terms of signal-to-noise ratio (SNR), rate of support recovery, Hamming error, and Hamming distance. In particular, OBCSA_CAF achieves 0.5–4 dB gains in SNR and [math]–[math] improvements in support recovery rate at low flipping ratios, whereas OBCSA_CAF_HOW yields 2–6.7 dB gains and [math]–[math] improvements, respectively, under high flipping ratios. Experiments on the MNIST and CIFAR-10 datasets further illustrate the superior performance of our methods in image reconstruction.

Recovering the Grating Profile from Limited-Aperture Observation: Data Retrieval and Shape Reconstruction

Tian Niu — 2026-05-08T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 1015-1036, June 2026.
Abstract.This paper proposes a two-stage framework to address the inverse diffraction grating problem with limited-aperture data. The first stage introduces a deep learning network for data retrieval, featuring a dual-branch, cross-attention architecture. Motivated by an information-theoretic analysis, this design is tailored to separate and adaptively fuse the diffracted field’s low- and high-frequency components, effectively handling their distinct noise sensitivities. The second stage employs a computationally efficient Newton-type algorithm for shape reconstruction, which avoids the need for a forward solver at each iteration. Numerical experiments show that our framework provides accurate and robust reconstructions.

Sampling Theory for Super-resolution with Implicit Neural Representations

Mahrokh Najaf — 2026-05-11T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 1037-1076, June 2026.
Abstract.Implicit neural representations (INRs) have emerged as a powerful tool for solving inverse problems in computer vision and computational imaging. INRs represent images as continuous domain functions realized by a neural network taking spatial coordinates as inputs. However, unlike traditional pixel representations, little is known about the sample complexity of estimating images using INRs in the context of linear inverse problems. Towards this end, we study the sampling requirements for recovery of a continuous-domain image from its low-pass Fourier samples by fitting a single hidden-layer INR with ReLU activation and a Fourier features layer using a generalized form of weight decay regularization. Our key insight is to relate minimizers of this nonconvex parameter space optimization problem to minimizers of a convex penalty defined over an infinite-dimensional space of measures. We identify a sufficient number of Fourier samples for which an image realized by an INR is exactly recoverable by solving the INR training problem. To validate our theory, we empirically assess the probability of achieving exact recovery of images realized by low-width single hidden-layer INRs and illustrate the performance of INRs on super-resolution recovery of continuous domain phantom images.

Deep Equilibrium Models for Poisson Imaging Inverse Problems via Mirror Descent

Christian Daniele — 2026-05-12T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 1077-1109, June 2026.
Abstract.Deep equilibrium models (DEQs) are implicit neural networks with fixed points that have recently gained attention for learning image regularization functionals, particularly in settings involving Gaussian fidelities, where assumptions on the forward operator ensure contractiveness of standard (proximal) gradient descent operators. In this work, we extend the application of DEQs to Poisson inverse problems, where the data fidelity term is more appropriately modeled by the Kullback–Leibler divergence. To this end, we introduce a novel DEQ formulation based on mirror descent defined in terms of a tailored non-Euclidean geometry that naturally adapts with the structure of the data term. This enables the learning of neural regularizers within a principled training framework. We derive sufficient conditions and establish refined convergence results based on the Kurdyka–Łojasiewicz framework for functions with nonclosed domains to guarantee the convergence of the learned reconstruction scheme and propose computational strategies that enable both efficient training and parameter-free inference. Numerical experiments show that our method outperforms traditional model-based approaches, and it is comparable to the performance of Bregman plug-and-play methods, while mitigating their typical drawbacks, such as time-consuming tuning of hyperparameters. The code is publicly available at https://github.com/christiandaniele/DEQ-MD.

Selective Focusing of Multiple Particles in a Layered Medium

Jun Lai — 2026-05-18T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 1110-1136, June 2026.
Abstract.Inverse scattering in layered media has a wide range of applications; examples include geophysical exploration, medical imaging, and remote sensing. In this paper, we develop a selective focusing method for identifying multiple unknown buried scatterers in a layered medium. The method is derived through an asymptotic analysis of the time reversal operator using the layered Green’s function and limited aperture measurements. We begin by showing the global focusing property of the time reversal operator. Then we demonstrate that each small sound-soft particle gives rise to one significant eigenvalue of the time reversal operator, while each sound-hard particle gives rise to three. The associated eigenfunction generates an incident wave that focuses selectively on the corresponding unknown particle. Finally, we employ the time reversal method as an initial indicator and propose an effective Bayesian inversion scheme to reconstruct multiple buried extended scatterers for enhanced resolution. Numerical experiments are provided to demonstrate the efficiency.

Topology-Guaranteed Image Segmentation: Enforcing Connectivity, Genus, and Width Constraints

Wenxiao Li — 2026-05-18T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 1137-1173, June 2026.
Abstract.Existing research highlights the crucial role of topological priors in image segmentation, particularly in preserving essential structures, such as connectivity and genus. Accurately capturing these topological features often requires incorporating width-related information, including the thickness and length inherent to the image structures. However, traditional mathematical definitions of topological structures lack this dimensional width information, limiting methods like persistent homology from fully addressing practical segmentation needs. To overcome this limitation, we propose a novel mathematical framework that explicitly integrates width information into the characterization of topological structures. This method leverages persistent homology, complemented by smoothing concepts from PDEs, to modify local extrema of upper level sets. This approach enables the resulting topological structures to inherently capture width properties. We incorporate this enhanced topological description into variational image segmentation models. Using some proper loss functions, we are also able to design neural networks that can segment images with the required topological and width properties. Through variational constraints on the relevant topological energies, our approach successfully preserves essential topological invariants, such as connectivity and genus counts, simultaneously ensuring that segmented structures retain critical width attributes, including line thickness and length. Numerical experiments demonstrate the effectiveness of our method, showcasing its capability to maintain topological fidelity while explicitly embedding width characteristics into segmented image structures.

Incomplete Data Multisource Static Computed Tomography Reconstruction with Diffusion Priors and Implicit Neural Representation

Ziju Shen — 2026-05-19T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 1174-1206, June 2026.
Abstract.The dose of X-ray radiation and the scanning time are crucial factors in computed tomography (CT) for clinical applications. In this work, we introduce a multisource static CT (MSCT) imaging system designed to rapidly acquire sparse view and limited angle data in CT imaging, addressing these critical factors. This linear imaging inverse problem is solved by a conditional generation process within the denoising diffusion image reconstruction framework. The noisy volume data sample generated by the reverse time diffusion process is projected onto the affine set to ensure its consistency with the measured data. To enhance the quality of the reconstruction, the 3D phantom’s orthogonal space projector is parameterized implicitly by a neural network. Then, a self-supervised learning algorithm is adopted to optimize the implicit neural representation. Through this multistage conditional generation process, we obtain a new approximate posterior sampling strategy for MSCT volume reconstruction. Numerical experiments are implemented with various imaging settings to verify the effectiveness of our methods for incomplete data MSCT volume reconstruction.

A Coarse-to-Fine Hybrid Registration and Fusion Framework for Hyperspectral Superresolution via Batch Image Alignment

Kunjing Yang — 2026-05-21T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 2, Page 1207-1243, June 2026.
Abstract.Fusing hyperspectral images (HSIs) with multispectral images (MSIs) has become a mainstream approach to enhance the spatial resolution of HSIs. Recently, many HSI-MSI fusion techniques have been developed and have achieved notable performance. Nevertheless, certain challenges still persist, including (a) Most existing fusion methods assume precise registration between HSIs and MSIs, which can be challenging to achieve in practice. (b) The obtained HSI-MSI pairs may not be fully utilized. To address these issues, we propose a coarse-to-fine hybrid registration and fusion approach. In the coarse stage, a joint registration and fusion (JRF) model is developed by integrating batch image alignment into the fusion process, which can achieve HSI-MSI registration and fusion simultaneously. To further enhance fusion performance, we propose a nonconvex low-rank and group-sparse (NLG) model in the fine stage. This model exploits the low-rank property to reconstruct the global image structure (i.e., the clean image), while simultaneously leveraging group sparsity to retrieve the textural details. This combination yields the proposed JRF-NLG method. Then, the JRF-NLG model is solved by the generalized Gauss–Newton (GGN) algorithm and the proximal alternating optimization (PAO) algorithm, respectively. Theoretically, we establish an error bound for the NLG model and prove the quadratic convergence rate of the GGN algorithm. Finally, extensive numerical experiments on simulated and real datasets are performed to verify the effectiveness of our method in registration and fusion. We also validate its effectiveness in enhancing classification performance.

A PCA Based Model for Surface Reconstruction from Incomplete Point Clouds

Hao Liu — 2026-01-02T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 1-34, March 2026.
Abstract.Point cloud data represents a crucial category of information for mathematical modeling, and surface reconstruction from such data is an important task across various disciplines. However, during the scanning process, the collected point cloud data may fail to cover the entire surface due to factors such as high light-absorption rate and occlusions, resulting in incomplete datasets. Inferring surface structures in data-missing regions and successfully reconstructing the surface poses a challenge. In this paper, we present a principal component analysis (PCA) based model for surface reconstruction from incomplete point cloud data. Initially, we employ PCA to estimate the normal information of the underlying surface from the available point cloud data. This estimated normal information serves as a regularizer in our model, guiding the reconstruction of the surface, particularly in areas with missing data. Additionally, we introduce an operator-splitting method to effectively solve the proposed model. Through systematic experimentation, we demonstrate that our model successfully infers surface structures in data-missing regions and well reconstructs the underlying surfaces, outperforming existing methodologies.

Diffusion at Absolute Zero: Langevin Sampling Using Successive Moreau Envelopes

Andreas Habring — 2026-01-02T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 35-77, March 2026.
Abstract.We propose a method for sampling from Gibbs distributions of the form [math] that leverages a family [math] of approximations of the target density which is deliberately constructed such that [math] exhibits favorable properties for sampling when [math] is large and such that [math] approaches [math] as [math] approaches 0. This sequence is obtained by replacing (parts of) the potential [math] with its Moreau envelope. Through the sequential sampling from [math] for decreasing values of [math] by a Langevin algorithm with appropriate step size, the samples are guided from a simple starting density to the more complex target quickly. We prove that [math] is Lipschitz continuous in the total variation distance and Hölder continuous in the Wasserstein-[math] distance, that the sampling algorithm is ergodic, and that it converges to the target density without assuming convexity or differentiability of the potential [math]. In addition to the theoretical analysis, we show experimental results that support the superiority of the method in terms of convergence speed and mode-coverage of multimodal densities over current algorithms. The experiments range from one-dimensional toy-problems to high-dimensional inverse imaging problems with learned potentials.

Analysis and Synthesis Denoisers for Forward-Backward Plug-and-Play Algorithms

Matthieu Kowalski — 2026-01-02T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 78-110, March 2026.
Abstract.In this work we study the behavior of the forward-backward (FB) algorithm when the proximity operator is replaced by a subiterative procedure to approximate a Gaussian denoiser, in a Plug-and-Play (PnP) fashion. Specifically, we consider both analysis and synthesis Gaussian denoisers within a dictionary framework, obtained by unrolling dual-FB iterations or FB iterations, respectively. We analyze the associated minimization problems as well as the asymptotic behavior of the resulting FB-PnP iterations. In particular, we show that the synthesis Gaussian denoising problem can be viewed as a proximity operator. For each case, analysis, and synthesis, we show that the FB-PnP algorithms solve the same problem whether we use only one or an infinite number of subiteration to solve the denoising problem at each iteration. To this aim, we show that each “one subiteration” strategy within the FB-PnP can be interpreted as a primal-dual algorithm when a warm-restart strategy is used. We further present similar results when using a Moreau–Yosida smoothing of the global problem, for an arbitrary number of subiterations. Finally, we provide numerical simulations to illustrate our theoretical results. In particular we first consider a toy compressive sensing example, as well as an image restoration problem in a deep dictionary framework.

Intrinsic Semiparametric Mixed Effects Model for Longitudinal Manifold-Valued Data

Di Xiong — 2026-01-02T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 111-139, March 2026.
Abstract.This paper introduces a novel intrinsic semiparametric mixed-effects model (ISMEM) to characterize the complex relationships between manifold-valued responses and Euclidean-valued covariates in the longitudinal setting. Manifold-valued data, characterized by their inherent nonlinearity, high dimensionality, and specific geometric structures, present significant challenges for longitudinal analysis. Most existing intrinsic models are developed for cross-sectional data and unsuited for longitudinal studies. To address this limitation, we propose ISMEM, an efficient framework designed for longitudinal manifold-valued datasets that capture such relationships at both group and individual levels. Key features of ISMEM include (i) the integration of fixed and random effects on Riemannian manifolds, enabling analysis at multiple levels; (ii) a semiparametric structure combining parametric and nonparametric components, enhancing both flexibility and interpretability; and (iii) the preservation of the geometric properties of manifold-valued observations, offering improved robustness and interpretability compared to Euclidean-based models. In addition, we develop a two-stage iterative estimation procedure and validate our approach through simulations on the symmetric positive definite (SPD) manifold. Finally, we apply ISMEM to longitudinal data from the Alzheimer’s Disease Neuroimaging Initiative (ADNI), reconstructing and comparing continuous three-dimensional (3D) shape trajectories of the lateral ventricle in Kendall shape space across distinct groups.

A Deep Neural Network Framework for Multivalued Mapping Problems with Varying Cardinality and Its Applications to Imaging

Geng Li — 2026-01-06T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 140-176, March 2026.
Abstract.This paper addresses the challenging problem of modeling and computing multivalued mappings with varying cardinality, where a single input can correspond to multiple valid outputs, and the number of possible outputs may vary for different inputs. Such scenarios are prevalent in various imaging applications where multiple plausible solutions exist for a given input. We introduce a deep neural network framework to model multivalued mappings with varying cardinality. The framework integrates a discrete codebook with a generative network to produce valid outputs for each input. The discrete codebook variables are combined with the input to guide the generator in producing different valid solutions. The discrete nature of the codebook enables the framework to efficiently estimate the conditional probability distribution of possible outputs through a fixed equiangular tight frame classifier. By jointly optimizing the discrete codebook and its uncertainty estimation during training using a specially designed covariance loss function, an accurate computation of multiple solution candidates with reliable confidence measures can be achieved. We demonstrate the effectiveness of the proposed framework on various imaging applications, using both synthetic and real datasets. Experimental results show the efficacy of our proposed model to generate multiple high-quality outputs while providing meaningful uncertainty estimates for each solution.

A Continuous Scale Space of Diffeomorphisms

Yechen Liu — 2026-01-06T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 177-206, March 2026.
Abstract.In this paper, we define and study a nested family of reproducing kernel Hilbert spaces of vector fields that is indexed by a range of scales, from which we construct a reproducing kernel Hilbert space of scale-dependent vector fields. We provide a characterization of the reproducing kernel of that space, with numerical approximations ensuring quick evaluations when this kernel does not have a closed form. We then introduce a multiscale version of the large deformation diffeomorphic metric mapping (LDDMM) problem and prove the existence of solutions. Finally, we provide numerical experiments performing landmark matching using multiscale LDDMM.

Spherical Area-Preserving Parameterization via Energy Minimization

Shu-Yung Liu — 2026-01-07T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 207-235, March 2026.
Abstract.We propose a novel method, called spherical authalic energy minimization (SAEM), for computing spherical area-preserving parameterizations of genus-zero closed surfaces with strong theoretical foundations. The global convergence of the associated computational algorithm is theoretically guaranteed. In addition, we introduce a Riemannian bijective correction method that ensures the bijectivity of the resulting mapping under mild assumptions. Numerical experiments show that SAEM effectively minimizes area distortion and achieves bijective mappings, outperforming state-of-the-art methods. Finally, we demonstrate the practical utility of SAEM in shape description.

Effective Solutions to Robust Orthogonal Nonnegative Matrix Factorization via Oblique Manifold Transformation

Fan Jia — 2026-01-08T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 236-265, March 2026.
Abstract.In this paper, we address the problem of robust orthogonal nonnegative matrix factorization (RONMF), a crucial challenge in data analysis. We first propose a RONMF model that explicitly handles both dense and sparse noise, making it suitable for a wide range of real-world applications. To circumvent the computational complexity associated with the Stiefel manifold and effectively solve the proposed model, we introduce an exact penalty method that transforms the optimization problem from the Stiefel manifold to the Oblique manifold. To achieve this, we develop the EP-RONMF algorithm, which seeks a point satisfying the weak second-order optimality conditions through an alternating proximal method and iterative updates of the penalty parameter. This approach successfully addresses the nonconvex nature of the ONMF problem and ensures convergence to a stable solution. To validate the efficacy of our method, we conducted extensive experiments on diverse datasets, including image, text, and hyperspectral data. The results clearly demonstrate the superiority of our approach compared to existing techniques.

Invertible ResNets for Inverse Imaging Problems: Competitive Performance with Provable Regularization Properties

Clemens Arndt — 2026-01-12T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 266-301, March 2026.
Abstract.Learning-based methods have demonstrated remarkable performance in solving inverse problems, particularly in image reconstruction tasks. Despite their success, these approaches often lack theoretical guarantees, which are crucial in sensitive applications such as medical imaging. Recent works by Arndt et al. addressed this gap by analyzing a data-driven reconstruction method based on invertible residual networks (iResNets). They revealed that, under reasonable assumptions, this approach constitutes a convergent regularization scheme. However, the performance of the reconstruction method was only validated on academic toy problems and small-scale iResNet architectures. In this work, we address this gap by evaluating the performance of iResNets on two real-world imaging tasks: a linear blurring operator and a nonlinear diffusion operator. To do so, we compare the performance of iResNets against state-of-the-art neural networks, revealing their competitiveness at the expense of longer training times. Moreover, we numerically demonstrate the advantages of the iResNet’s inherent stability and invertibility by showcasing increased robustness across various scenarios as well as interpretability of the learned operator, thereby reducing the black-box nature of the reconstruction scheme.

A Neural Network–Enhanced Born Approximation for Inverse Scattering

Ansh Desai — 2026-02-13T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 302-326, March 2026.
Abstract.Time-harmonic acoustic inverse scattering concerns the ill-posed and nonlinear problem of determining the refractive index of an inaccessible, penetrable scatterer based on far-field wave scattering data. When the scattering is weak, the regularized inverse Born approximation provides a linearized model for recovering the shape and material properties of a scatterer. We propose two convolutional neural network (CNN) algorithms to correct the traditional inverse Born approximation even when the scattering is not weak. These are denoted Born-CNN (BCNN) and CNN-Born (CNNB). BCNN applies a postcorrection to the Born reconstruction, while CNNB precorrects the data. Both methods leverage the Born approximation’s excellent fidelity in weak scattering while extending its applicability beyond its theoretical limits. CNNB particularly exhibits a strong generalization to more complex out-of-distribution scatterers. Based on numerical tests and benchmarking against other standard approaches, our corrected Born models provide alternative data-driven methods for obtaining the refractive index, extending the utility of the Born approximation to regimes where the traditional method fails.

Compression of Currents and Varifolds

Allen Paul — 2026-02-25T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 327-363, March 2026.
Abstract.We derive an algorithm for compression of the currents and varifolds representations of shapes using Ridge Leverage Score sampling and the theory of Nystrom approximation in Reproducing Kernel Hilbert Spaces. Our method is faster than existing compression techniques and comes with theoretical guarantees on the rate of convergence of the compressed approximation as a function of the smoothness of the associated shape representation. The obtained compressions are shown to be useful for down-line tasks such as nonlinear shape registration in the Large Deformation Diffeomorphic Metric Mapping (LDDMM) framework, even for very high compression ratios. The performance of our algorithm is demonstrated on large-scale shape data from modern geometry processing datasets and is shown to be fast and scalable with rapid error decay.

Subspace Method of Moments for Ab Initio 3-D Single Particle Cryo-EM Reconstruction

Jeremy Hoskins — 2026-02-27T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 364-409, March 2026.
Abstract.Cryo-electron microscopy (cryo-EM) is a widely used technique for recovering the three-dimensional (3-D) structure of biological molecules from a large number of experimentally generated noisy 2-D tomographic projection images of the 3-D structure, taken from unknown viewing angles. Through computationally intensive algorithms, these observed images are processed to reconstruct the 3-D structures. Many popular computational methods rely on estimating the unknown angles as part of the reconstruction process, which becomes particularly challenging at low signal-to-noise ratios. The method of moments offers an alternative approach that circumvents the estimation of viewing orientations of individual projection images by instead estimating the underlying distribution of the viewing angles, and is robust to noise given sufficiently many images. However, the method of moments typically entails computing higher-order moments of the projection images, incurring significant computational and memory costs. To mitigate this, we propose a new approach called the subspace method of moments (SubspaceMoM), which compresses the first three moments using data-driven low-rank tensor techniques as well as expansion into a suitable function basis. The compressed moments can be efficiently computed from the set of projection images using numerical quadrature and can be employed to jointly reconstruct the 3-D structure and the distribution of viewing orientations. We illustrate the practical applicability of SubspaceMoM through numerical experiments using up to the third-order moment on synthetic datasets with a simplified cryo-EM image formation model, which significantly improves the reconstruction resolution compared to previous MoM approaches.

Refining Image Edge Detection via Linear Canonical Riesz Transforms

Shuhui Yang — 2026-03-04T08:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 410-445, March 2026.
Abstract.By combining the linear canonical transform and the Riesz transform, we introduce the linear canonical Riesz transform (LCRT), which is further proved to be a linear canonical multiplier. Using this LCRT multiplier, we conduct numerical simulations on images. Notably, due to the fact that the linear canonical transform itself admits a fast algorithm, the LCRT can achieve a computational speed comparable to that of the linear canonical transform. Based on this, we introduce the new concept of the sharpness [math] of the edge strength and continuity of images associated with the LCRT and, using it, we propose a new LCRT image edge detection method (LCRT-IED method) and provide its mathematical foundation. Our experiments indicate that this sharpness [math] characterizes the macroscopic trend of edge variations of the image under consideration, while this new LCRT-IED method not only controls the overall edge strength and continuity of the image, but also excels in feature extraction in some local regions. These highlight the fundamental differences between the LCRT and the Riesz transform, which are precisely due to the multiparameter of the former. This new LCRT-IED method might be of significant importance for image feature extraction, image matching, and image refinement.

Accurate, Provable, and Fast Polychromatic Tomographic Reconstruction: A Variational Inequality Approach

Mengqi Lou — 2026-03-09T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 446-479, March 2026.
Abstract.We consider the problem of signal reconstruction for computed tomography (CT) given a nonlinear forward model that accounts for exponential signal attenuation, a polychromatic X-ray source, general measurement noise (e.g., Poisson shot noise), and observations acquired over multiple wavelength windows. We develop a simple iterative algorithm for single-material reconstruction, which we call EXACT (EXtragradient Algorithm for Computed Tomography), based on formulating our estimate as the fixed point of a monotone variational inequality. We prove guarantees on the statistical and computational performance of EXACT given realistic assumptions on the measurement process. We also consider a recently introduced variant of this model with Gaussian measurements and present sample and iteration complexity bounds for EXACT that improve upon those of existing algorithms. We apply our EXACT algorithm to a CT phantom image recovery task and show that it often requires fewer X-ray views, lower source intensity, and less computation time to achieve similar reconstruction quality to existing methods. Code is available at https://github.com/voilalab/exact.

Robust Acoustic and Elastic Full Waveform Inversion by Adaptive Tikhonov-TV Regularization

Kamal Aghazade — 2026-03-10T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 480-518, March 2026.
Abstract.Full waveform inversion (FWI) is a powerful wave-based imaging technique, but its inherent ill-posedness and nonconvexity make it prone to local minima and poor convergence. Regularization techniques are commonly employed to stabilize FWI by incorporating prior information that enforces structural constraints, such as smooth variations or piecewise-constant behavior. Among them, Tikhonov regularization promotes smoothness, while total variation (TV) regularization preserves sharp boundaries; both are widely used in solving ill-posed inverse problems. However, in the context of FWI, we highlight two key shortcomings of these regularization methods. First, subsurface model parameters (P- and S-wave velocities, density) often exhibit complex geological formations with sharp discontinuities separating distinct layers, while parameters within each layer vary smoothly. Neither Tikhonov nor TV regularization alone can effectively constrain such piecewise-smooth structures. Second, and more critically, when the initial model is far from the true model, these regularization assumptions can lead to convergence toward a local minimum. To address these limitations, we propose an adaptive Tikhonov-TV (TT) regularization method that decomposes the model into smooth and blocky components, enabling robust recovery of piecewise-smooth structures. The method is implemented within the alternating direction method of multipliers (ADMM) framework and incorporates an automated balancing strategy based on robust statistical analysis. Extensive numerical experiments on both acoustic and elastic FWI are conducted using challenging benchmark geological models. The results demonstrate that TT regularization significantly improves convergence and reconstruction accuracy compared to Tikhonov and TV regularization when applied separately. We show that for complex models and remote initial models, both Tikhonov and TV regularization tend to converge to local minima, whereas the TT regularization effectively mitigates cycle skipping through its adaptive combination of the two regularization strategies.

Shape Prior Segmentation Guided by Harmonic Beltrami Signature

Chenran Lin — 2026-03-12T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 519-554, March 2026.
Abstract.This paper presents a novel shape prior segmentation model guided by the harmonic Beltrami signature (HBS) that integrates shape information to enhance the segmentation accuracy. The HBS is a robust shape representation that fully captures a 2D simply connected shape. It exhibits resilience to perturbations and is invariant under translation, rotation, and scaling. These properties make it a suitable candidate for incorporating general shape information into the segmentation model, instead of relying on primitive properties such as convexity or topology. Our shape prior segmentation model embeds the HBS within a quasi-conformal, topology-preserving segmentation framework. With the shape prior knowledge, our proposed model significantly enhances the segmentation performance, especially for low-quality or occluded images. The robustness of HBS to shape perturbations allows for the use of a simple [math] distance metric to define shape dissimilarity, which simplifies the optimization process. Besides, the invariance of HBS to translation, rotation, and scaling allows the integration of shape prior information without requiring the knowledge of orientation, size and position of the shape in the image. This enhances its practicality for real-world applications. Extensive experiments have been carried out on both synthetic and real images. Results demonstrate that our proposed model significantly improves segmentation accuracy, particularly for degraded or corrupted images.

Quasi-conformal Convolution: A Learnable Convolution for Deep Learning on Simply Connected Open Surfaces

Han Zhang — 2026-03-13T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 555-583, March 2026.
Abstract.Deep learning on non-Euclidean domains is important for analyzing complex geometric data that lacks common coordinate systems and familiar Euclidean properties. A central challenge in this field is to define convolution on domains, which inherently possess irregular and non-Euclidean structures. In this work, we introduce quasi-conformal convolution (QCC), a novel framework for defining convolution on simply-connected open surfaces using quasi-conformal theories. Each QCC operator is linked to a specific quasi-conformal mapping, enabling the adjustment of the convolution operation through manipulation of this mapping. By utilizing trainable estimator modules that produce quasi-conformal mappings, QCC facilitates adaptive and learnable convolution operators that can be dynamically adjusted according to the underlying data structured on the surfaces. QCC unifies a broad range of spatially defined convolutions, facilitating the learning of tailored convolution operators on each underlying surface optimized for specific tasks. Building on this foundation, we develop the quasi-conformal convolutional neural network (QCCNN) to address a variety of tasks related to geometric data. We validate the efficacy of QCCNN through the classification of images defined on curvilinear simply-connected open Riemann surfaces, demonstrating superior performance in this context. Additionally, we explore its potential in medical applications, including craniofacial analysis using 3D facial data and lesion segmentation on 3D human faces, achieving enhanced accuracy and reliability.

Total Generalized Variation of the Normal Vector Field and Applications to Mesh Denoising

Lukas Baumgärtner — 2026-03-17T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 584-611, March 2026.
Abstract.We propose a novel formulation for the second-order total generalized variation (TGV) of the normal vector on an oriented, triangular mesh embedded in [math]. The normal vector is considered as a manifold-valued function, taking values on the unit sphere. Our formulation extends previous discrete TGV models for piecewise constant scalar data that utilize a Raviart–Thomas function space. To extend this formulation to the manifold setting, a tailor-made tangential Raviart–Thomas-type finite element space is constructed in this work. The new regularizer is compared to existing methods in mesh denoising experiments.

An Inverse Obstacle Scattering Problem with Passive Data in the Time Domain

Xiaoli Liu — 2026-03-19T07:00:00Z

SIAM Journal on Imaging Sciences, Volume 19, Issue 1, Page 612-642, March 2026.
Abstract.This work considers a time domain inverse acoustic obstacle scattering problem due to passive data. Motivated by the Helmholtz–Kirchhoff identity in the frequency domain, we propose to relate the time domain measurement data in passive imaging to an approximate data set given by the subtraction of two scattered wave fields. We propose a time domain linear sampling method for the approximate data set and show how to tackle the measurement data in passive imaging. An imaging functional is built based on the linear sampling method, which reconstructs the support of the unknown scattering object using directly the time domain measurements. The functional framework is based on the Laplace transform, which relates the mapping properties of Laplace domain factorized operators to their counterparts in the time domain. Numerical examples are provided to illustrate the capability of the proposed method.