Modified Papoulis-Gerchberg algorithm for sparse signal recovery

Motivated by the well-known Papoulis-Gerchberg algorithm, an iterative thresholding algorithm for recovery of sparse signals from few observations is proposed. The sequence of iterates turns out to be similar to that of the thresholded Landweber iterations, although not the same. The performance of the proposed algorithm is experimentally evaluated and compared to other state-of-the-art methods.

💡 Research Summary

The paper introduces a novel iterative thresholding algorithm for recovering sparse signals from a severely undersampled set of linear measurements. Inspired by the classic Papoulis‑Gerchberg method, which alternates between enforcing a known‑segment constraint and a band‑limiting constraint for continuous‑time signals, the authors adapt the core idea to the discrete, sparsity‑constrained setting. Their algorithm consists of two distinct operations per iteration: (1) a sparsity‑enforcing soft‑thresholding step applied to the current signal estimate, and (2) a measurement‑consistency step that projects the thresholded estimate back onto the space defined by the linear measurement operator using its adjoint. Formally, with K denoting the measurement matrix and K* its adjoint, the iteration can be written as
 fₙ = K* ( g + Sγ(fₙ₋₁) – K Sγ(fₙ₋₁) ), f₀ = K* g,
where Sγ is the soft‑thresholding operator with threshold γ. This formulation resembles the well‑known thresholded Landweber iterations but differs in that the thresholding is performed before the back‑projection, explicitly separating “constraint maintenance” from “observation recovery.”

The authors do not provide a formal convergence proof; instead they rely on existing results for soft‑thresholded Landweber schemes, suggesting that similar convergence behavior is plausible. The algorithm’s simplicity is a notable advantage: it requires only matrix‑vector multiplications with K and K* and a pointwise thresholding operation, and the threshold γ can be automatically selected using the Birge‑Massart rule applied to stationary wavelet transform (SWT) coefficients (Haar wavelet, one level).

Experimental validation is carried out on both synthetic 1‑D signals and 2‑D medical‑imaging data. For the 1‑D HeaviSine benchmark (N = 1024), the method is tested with 70, 100, 150, and 200 randomly selected samples. Mean‑squared error (MSE) results show consistent improvement over two state‑of‑the‑art baselines: ℓ₁‑norm minimization (basis pursuit) and total‑variation (TV) minimization. For example, with 150 samples the proposed method achieves an MSE of 6.5 × 10⁻³, compared to 3.8 × 10⁻² (ℓ₁) and 1.31 × 10⁻² (TV).

In the 2‑D experiments, a 256 × 256 Shepp‑Logan phantom is reconstructed from highly undersampled Fourier data obtained along K = 9, 11, 15, 21 radial lines. Reconstruction quality is measured by peak signal‑to‑noise ratio (PSNR). The proposed algorithm yields PSNR values of 24.97 dB (K = 9), 39.31 dB (K = 15), and 199.75 dB (K = 21), substantially higher than the ℓ₁ baseline (14.36 dB – 21.13 dB) and the TV baseline (11.8 dB – 27.1 dB). Visual inspection confirms that the proposed reconstructions retain fine structures and exhibit far fewer artifacts.

A further set of experiments adds white Gaussian noise to the measurements. Even under significant noise levels, the algorithm produces reconstructions whose PSNR exceeds that of the noisy measurements themselves, demonstrating robustness to measurement perturbations.

Overall, the paper makes three key contributions: (1) a clear adaptation of the Papoulis‑Gerchberg framework to sparse signal recovery, (2) a simple yet effective iterative scheme that blends soft‑thresholding with a Landweber‑type back‑projection, and (3) empirical evidence that the method outperforms leading ℓ₁ and TV approaches on both 1‑D and 2‑D tasks, including noisy scenarios. The authors acknowledge that a rigorous convergence analysis and extensions to hard‑thresholding or alternative transforms are left for future work, suggesting a promising research direction for the community.