Convex Data-Fitting Formulation

• Convex data-fitting formulation is an optimization technique that minimizes a strictly convex objective by integrating quadratic fidelity with non-smooth regularizers.
• It employs quadratic data terms and regularization methods like vector total variation to enforce measurement consistency and preserve key image features.
• The approach utilizes variable splitting and ADMM to efficiently solve large-scale inverse problems while ensuring robust convergence and computational scalability.

A convex data-fitting formulation refers to an optimization approach in which the parameters of a model (often a high-dimensional signal or function) are estimated by minimizing a strictly convex objective function. This objective incorporates quadratic data-fitting terms, modeling measurement consistency with observed data under linear degradations and noise, as well as regularization terms—possibly non-quadratic and non-smooth, such as vectorial total variation (VTV)—to encode desirable structure in the solution. Convexity ensures that any local minimum is a global minimum, enabling the use of efficient and reliable optimization algorithms. This paradigm is particularly impactful in ill-posed inverse problems, where direct inversion is unstable or intractable due to conditioning, dimensionality, or structural complexity.

1. Mathematical Structure of Convex Data-Fitting Formulations

Convex data-fitting problems, as exemplified by hyperspectral image superresolution, are constructed by combining quadratic data fidelity terms with regularization:

$$\min_X \frac{1}{2} \|Y_h - E X B M\|_F^2 + \frac{\lambda_m}{2} \|Y_m - R E X\|_F^2 + \lambda_\phi \phi(X)$$

• $Y_h$ denotes the observed hyperspectral data (high spectral, low spatial resolution).
• $Y_m$ represents the multispectral data (low spectral, high spatial resolution).
• $E \in \mathbb{R}^{L_h \times s}$ is a basis or subspace matrix capturing the low-dimensional structure of hyperspectral data; the high-resolution image is parameterized as $Z = EX$.
• $B$ and $M$ are the spatial blurring and downsampling operators, respectively.
• $R$ is the spectral response matrix.
• $\phi(X)$ is a regularization term imposing spatial smoothness and edge alignment among hyperspectral bands, often realized as vector total variation (VTV):

$$\phi(X) = \sum_{j=1}^{n_m} \sqrt{ \sum_{i=1}^{s} \left[ (X D_h)_{ij}^2 + (X D_v)_{ij}^2 \right] }$$

where $D_h$ and $D_v$ are discrete horizontal and vertical finite-difference operators.

The quadratic terms model the fidelity to observations under the degradations and noise introduced by the sensing process, while the regularizer encodes prior information about the structure of the solution, such as spatial coherence and edge preservation.

2. Variable Splitting and ADMM for Optimization

The convex objective is generally challenging to optimize directly due to the presence of non-diagonalizable operators (e.g., downsampling $M$) and non-smooth regularization (e.g., VTV). To address this, the problem is reformulated using variable splitting:

$$\begin{aligned} &\min_{X, V_1, V_2, V_3, V_4} \frac{1}{2} \| Y_h - E V_1 M\|_F^2 + \frac{\lambda_m}{2} \|Y_m - R E V_2\|_F^2 + \lambda_\phi \phi(V_3, V_4) \ &\text{subject to} \quad X B = V_1, \quad X = V_2, \quad X D_h = V_3, \quad X D_v = V_4 \end{aligned}$$

This decouples the strongly-coupled terms, allowing the use of specialized proximal or linear solvers for each block. The Alternating Direction Method of Multipliers (ADMM) is then applied, alternating minimization over the variables $X, V_1, ..., V_4$, and dual variable updates:

• Updates for $X$, $V_1$, and $V_2$ are quadratic or involve (block-)circulant operators, making them tractable using FFT-based inversion.
• Updates for $V_3$ and $V_4$ are solved efficiently via the proximal mapping of the VTV norm—implemented as pixel-wise vector soft-thresholding (proximal operator of the $\ell_{2,1}$ norm).

ADMM robustly converges under mild conditions, and the variable splitting isolates the computational and structural bottlenecks, making large-scale problems tractable.

3. Modeling Degradations and Subspace Structure

A critical modeling aspect is the representation $Z = E X$, where $E$ spans a subspace capturing the spectral correlation structure of hyperspectral images. This subspace, often determined by techniques such as vertex component analysis (VCA) or truncated SVD, dramatically reduces dimensionality and estimation complexity (typical reduction: $L_h \sim 200$ spectral bands to $s \sim 10$ subspace components).

The observation models are:

$$\begin{aligned} Y_h &= Z B M + N_h \ Y_m &= R Z + N_m \end{aligned}$$

where $N_h, N_m$ are additive noise.

By projecting physical measurements into a low-dimensional subspace and formulating the fusion in the coordinate space $X$, the method leverages prior knowledge about the signal structure and increases both statistical efficiency and numerical stability.

4. Handling Non-Diagonalizable Operators and Nonsmooth Regularization

The inclusion of operators such as $M$ (downsampling) impedes direct inversion by Fourier techniques since such matrices do not commute with the shift operators. The introduced variable splits are chosen explicitly to localize these non-invertible operators' effects; FFT-based solvers can be exploited where possible (for $B$ and $D_h, D_v$ with circular boundary conditions).

The nonsmooth, non-quadratic nature of the VTV regularizer is overcome by proximal splitting: the update for $V_3, V_4$ is executed using the closed-form vector soft thresholding operator, which, for each spatial location, shrinks the joint gradient vector towards zero up to a radius determined by $\lambda_\phi$. This operation strongly promotes spatial smoothness while preserving sharp transitions—i.e., aligned image edges—across bands.

5. Numerical Performance and Scalability

The convex data-fitting approach, together with subspace modeling, yields highly scalable algorithms. Experimental results indicate:

• State-of-the-art accuracy in superresolution benchmarks: lower ERGAS values, reduced Spectral Angle Mapper (SAM), and higher Universal Image Quality Index (UIQI) compared to pansharpening and other fusion methods. These metrics directly quantify reconstruction error, spectral fidelity, and perceptual visual quality.
• Robust convergence typically within 200 iterations of ADMM/SALSA.
• Subspace reduction from hundreds to $\sim$ tens of variables per spatial pixel enables efficient handling of extremely large images.
• Unlike other methods that require diagonalizable degradations or regularizers, this formulation is robust to arbitrary downsampling and non-differentiable structure in the regularization.

The computational strategy—splitting, ADMM, and FFT acceleration—enables tractable inference for massive-scale imaging data, as frequently encountered in remote sensing.

6. Broader Applicability and Extensions

The convex data-fitting formulation with structured regularization and variable splitting is broadly applicable beyond hyperspectral fusion:

• It has been adopted in image deblurring, denoising, medical imaging (e.g., MRI and CT reconstruction), compressed sensing, and other inverse problems requiring global optimality and structural constraints.
• The design principles allow adaptation to domain-specific degradations, arbitrary convex but possibly nonsmooth regularizers, and diverse structural priors (anisotropic TV, nuclear norm, group sparsity).
• The decoupling and proximal splitting techniques are extensible to more general ADMM or primal-dual proximal frameworks, facilitating integration with learned components or plug-and-play priors, provided convexity is preserved within the core formulation.

The rigorous separation of fidelity and regularization, explicit modeling of all degradations, and the direct enabling of principled optimization with strong guarantees underpin the distinctive power of the convex data-fitting approach in large-scale, ill-posed, and structured signal estimation problems (Simões et al., 2014).