Monte Carlo Maximum Likelihood Reconstruction for Digital Holography with Speckle

In coherent imaging, speckle is statistically modeled as multiplicative noise, posing a fundamental challenge for image reconstruction. While maximum likelihood estimation (MLE) provides a principled framework for speckle mitigation, its application to coherent imaging system such as digital holography with finite apertures is hindered by the prohibitive cost of high-dimensional matrix inversion, especially at high resolutions. This computational burden has prevented the use of MLE-based reconstruction with physically accurate aperture modeling. In this work, we propose a randomized linear algebra approach that enables scalable MLE optimization without explicit matrix inversions in gradient computation. By exploiting the structural properties of sensing matrix and using conjugate gradient for likelihood gradient evaluation, the proposed algorithm supports accurate aperture modeling without the simplifying assumptions commonly imposed for tractability. We term the resulting method projected gradient descent with Monte Carlo estimation (PGD-MC). The proposed PGD-MC framework (i) demonstrates robustness to diverse and physically accurate aperture models, (ii) achieves substantial improvements in reconstruction quality and computational efficiency, and (iii) scales effectively to high-resolution digital holography. Extensive experiments incorporating three representative denoisers as regularization show that PGD-MC provides a flexible and effective MLE-based reconstruction framework for digital holography with finite apertures, consistently outperforming prior Plug-and-Play model-based iterative reconstruction methods in both accuracy and speed. Our code is available at: https://github.com/Computational-Imaging-RU/MC_Maximum_Likelihood_Digital_Holography_Speckle.

💡 Research Summary

This paper addresses the long‑standing challenge of speckle noise in coherent imaging systems, particularly digital holography, by formulating a physically accurate maximum‑likelihood estimation (MLE) problem that explicitly models speckle as multiplicative noise. The forward model is y = A g + z, where g ∼ 𝒞𝒩(0, diag(x)) represents the complex reflectance with speckle, x is the underlying speckle‑free reflectivity, A = F⁻¹ M F encodes the optical aperture (circular or annular) in the Fourier domain, and z is additive Gaussian sensor noise. The negative log‑likelihood f(x) = log det Σ(x) + yᴴΣ⁻¹(x)y, with Σ(x) = A diag(x)Aᴴ + σ²_z I, is convex in x but requires the inversion of a large m × m covariance matrix, which is computationally prohibitive for high‑resolution images.

The authors propose a novel algorithm called Projected Gradient Descent with Monte‑Carlo estimation (PGD‑MC) that eliminates the need for explicit matrix inversions while preserving the exact likelihood. Two key ideas enable this: (1) The term AᴴΣ⁻¹(x) y is computed by solving Σ(x) h = y with the conjugate‑gradient (CG) method, exploiting the Hermitian positive‑definite nature of Σ(x). (2) The diagonal of AᴴΣ⁻¹(x)A, required for the gradient, is estimated via Monte‑Carlo trace estimation. Random Gaussian vectors v are drawn, the linear systems Σ(x) u = A v are solved with CG, and w = Aᴴu is formed; the element‑wise product w⊙v averaged over K samples yields an unbiased estimator of the diagonal. This stochastic approximation avoids any explicit formation of Σ⁻¹ or AᴴΣ⁻¹A.

With the gradient ∇f_L(x_t) = d_t − (1/L)∑ℓ|m_ℓ|² (where d_t is the Monte‑Carlo diagonal estimate and m_ℓ = Aᴴh_ℓ), the algorithm performs a standard projected gradient step: x{t+1} = Π_X(x_t − μ∇f_L(x_t)). The projection Π_X encodes prior knowledge through image denoisers; the authors demonstrate three choices—BM3D, DnCNN, and Deep Image Prior (DIP)—showing that the framework is agnostic to the specific regularizer. Because A is a composition of FFT, IFFT, and a binary mask, all matrix‑vector products are implemented with O(n log n) complexity using fast Fourier transforms, allowing the method to scale to 512 × 512 images without excessive memory consumption.

Experimental results cover both circular and annular apertures, three speckle variance levels (σ_z = 1.5, 2.5, 5), and varying numbers of measurements L = 4. Across all settings, PGD‑MC consistently outperforms state‑of‑the‑art Plug‑and‑Play (PnP) methods such as CP‑nP‑EM and traditional EM‑based approaches. Quantitatively, PGD‑MC achieves 2–4 dB higher PSNR and reduces reconstruction time by a factor of three or more. The Monte‑Carlo sample size K = 5–10 and CG iteration limit of 10–20 are sufficient for stable convergence, demonstrating the practical efficiency of the stochastic diagonal estimator.

The main contributions are: (i) a mathematically exact MLE formulation that retains realistic aperture modeling; (ii) a randomized linear‑algebra plus CG scheme that bypasses high‑dimensional matrix inversions; (iii) a flexible PGD‑MC framework that can incorporate any denoiser‑based prior; and (iv) extensive empirical validation showing superior accuracy and speed over existing model‑based reconstructions. The authors suggest future extensions to non‑linear or multi‑layer aperture models, multi‑wavelength holography, and hardware acceleration (GPU/FPGA) for real‑time deployment. Overall, PGD‑MC represents a significant step toward practical, high‑quality speckle mitigation in digital holography and other coherent imaging modalities.