Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Maximum-likelihood estimation in ptychography in the presence of Poisson–Gaussian noise statistics

Open Access Open Access

Abstract

Optical measurements often exhibit mixed Poisson–Gaussian noise statistics, which hampers the image quality, particularly under low signal-to-noise ratio (SNR) conditions. Computational imaging falls short in such situations when solely Poissonian noise statistics are assumed. In response to this challenge, we define a loss function that explicitly incorporates this mixed noise nature. By using a maximum-likelihood estimation, we devise a practical method to account for a camera readout noise in gradient-based ptychography optimization. Our results, based on both experimental and numerical data, demonstrate that this approach outperforms the conventional one, enabling enhanced image reconstruction quality under challenging noise conditions through a straightforward methodological adjustment.

© 2023 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement

Corrections

16 November 2023: A typographical correction was made to the Fig. 1 caption.

In the rapidly evolving field of computational imaging, ptychography has emerged as a powerful technique capable of producing high-resolution phase and amplitude images from diffraction patterns. It involves translating a thin object through overlapping illuminations and measuring the resulting diffraction patterns behind the object with a camera sensor [1,2]. Subsequently, the complex-valued image is constructed through an iterative optimization algorithm, necessitating the formulation and minimization of a loss function, alternatively referred to as the objective, cost, or error function. Ptychography has found applications in a wide range of topics, including label-free biological imaging [35], optical metrology [610], and atomic-resolution imaging using electron beams [1113].

Fundamentally, the basis of ptychographic reconstructions is the detection of photon counts on a camera sensor and, therefore, subject to Poissonian noise even under ideal measurement conditions. Given the assumption of an underlying noise model, a powerful and robust optimization strategy is the maximum-likelihood estimation (MLE) principle [14]. By leveraging MLE in ptychography, one seeks to estimate the studied object parameters that render the observed diffraction patterns most probable [1519]. However, an additive camera readout noise is often neglected, thus leaving a gap in the fidelity of reconstructions. This is a salient concern as the presence of readout noise, typically Gaussian, is an important element of practical ptychographic measurements within the visible spectrum. Ignoring this noise source oversimplifies the underlying statistical model and introduces errors to the reconstructed image, especially when the detected photon counts and the signal-to-noise ratio (SNR) are low. The distinction between a Poissonian and a mixed Poisson–Gaussian noise model is depicted in the simulated images presented in Fig. 1.

 figure: Fig. 1.

Fig. 1. Comparative visualization of a grayscale image distorted by different noise types. Panel (A) depicts an image with simulated Poissonian noise, while panel (B) illustrates the effects of simulated mixed Poisson–Gaussian noise resulting from additive readout noise. Inlaid values indicate peak signal-to-noise ratio (PSNR) with respect to the ground truth.

Download Full Size | PDF

In this Letter, we propose a loss function for automatic differentiation ptychography that explicitly incorporates both Poissonian and Gaussian noise sources. This approach brings us closer to the real-world conditions of ptychographic measurements, thereby paving the way for superior performance in image reconstruction under challenging noise conditions. We outline a practical method to incorporate camera readout noise in computational imaging. Furthermore, we provide a comprehensive comparison between the image reconstruction quality using a mixed-statistics loss function and that of a conventional loss function which presumes solely Poissonian noise statistics. For this, we present reconstruction results obtained from both experimental and numerical data.

In ptychography, the typical reconstruction approach involves minimizing a loss function representing the difference between the intensity values of the observed diffraction pattern $X_k$ and the anticipated diffraction patterns $I_k(\boldsymbol{\theta })$ as determined by a parameter set $\boldsymbol{\theta }$, which embodies the object under investigation, at all pixel locations indexed by $k$. In the presence of measurement noise, it is insightful to tackle the problem of ptychographic reconstruction by maximizing the likelihood of the observed given the object parameters. From this probabilistic perspective, one seeks the object parameters that make the observed data most likely, which renders object retrieval more robust in scenarios of low SNR.

In the Supplement 1 (section 1), we elaborate on deriving the maximum-likelihood estimation (MLE) loss functions for two different types of noise statistics. When operating under the assumption of Poissonian counting noise, the loss function $L_{\mathrm {Poisson}}$ that yields the maximum-likelihood estimate is expressed as

$$L_{\mathrm{Poisson}}(\boldsymbol{\theta})=\sum_{k=1}^N\left( \sqrt{X_k}- \sqrt{I_k(\boldsymbol{\theta})}\right)^2,$$
where the sum encompasses $N$ statistically independent pixels on the camera sensor. To account for Gaussian readout noise on the camera sensor, an additional data acquisition step becomes essential to extend the MLE loss function: the variance $\sigma _k^2$ of the readout noise of the camera must be determined using multiple full-frame dark images. With this additional information and the assumption that the Poissonian component of the statistics can be approximated by a Gaussian distribution, we formulate the MLE loss function that incorporates mixed Poisson–Gaussian noise statistics:
$$L_{\mathrm{Mixed}}(\boldsymbol{\theta})=\sum_{k=1}^N\left( \ln[I_k(\boldsymbol{\theta}) + \sigma_k^2] + \frac{[X_k- I_k(\boldsymbol{\theta})]^2}{I_k(\boldsymbol{\theta}) + \sigma_k^2} \right).$$
Note that this expression is not only relevant for pixels with high photon counts where the Gaussian approximation of Poisson statistics is most accurate but also for pixels with low photon counts. Indeed, for those with low-count pixels, the Gaussian readout noise is the dominant source of noise such that the deviation of Poisson statistics from a Gaussian distribution becomes irrelevant.

To validate the beneficial effect of the loss function $L_{\mathrm {Mixed}}(\boldsymbol{\theta })$ experimentally, we are considering a standard ptychography setup in a transmission geometry (Fig. 2). A circular 500-µm pinhole is illuminated with coherent light with a wavelength of 561 nm and relayed to the sample plane using two lenses with a magnification of $M=3$. Therefore, a binary target sample is illuminated at 80 scanning positions with an overlap of approximately 60 % between adjacent positions, which lies within the ideal regime according to [20]. The scattered light is captured by a CMOS camera positioned 38 mm away from the object. To ensure a comprehensive comparison as a function of SNR, we employ four different exposure settings per scanning position, spanning a range from 30 µs to 300 ms, with each subsequent exposure time differing by a factor of 10 (see Fig. 3). For each exposure time, we determine the variances $\sigma _k^2$ by capturing a stack of 300 dark images. This quantifies the readout noise level associated with each pixel $k$ of our camera. A comprehensive overview of further details about the experimental implementation and methodology can be found in Section 2 of Supplement 1.

 figure: Fig. 2.

Fig. 2. Schematic drawing of our ptychography setup used in the experiment and for numerical simulations. A 500-µm pinhole is illuminated and relayed onto the object using a 2-lens system. The object is shifted laterally through the beam using an XY stage. A CMOS camera records the diffraction intensities at a distance of 38 mm downstream of the object.

Download Full Size | PDF

 figure: Fig. 3.

Fig. 3. Comparison of image reconstruction qualities from ptychographic datasets with decreasing signal-to-noise ratio (SNR) from left to right. Top row: Ptychographic datasets at different camera exposures with total illumination intensity, respectively. Second row: Reconstructions using a loss function assuming solely Poissonian noise statistics. Third row: Reconstructions using the loss function defined by Eq. (2), which incorporates mixed Poisson–Gaussian noise statistics. The edge length of every image equals 3.5 mm.

Download Full Size | PDF

By maintaining a constant illumination power, we acquire four distinct ptychographic datasets, each corresponding to a different signal-to-noise ratio (SNR), as demonstrated in the top row of Fig. 3. With these datasets in hand, we proceed to perform image reconstructions utilizing two different loss functions (Eqs. (1) and (2)) within an automatic differentiation-based ptychography framework, as detailed in [21] and similar to [22]. For this analysis, we precalibrate the illumination field using an additional high SNR measurement and restrict our optimizations to the complex-valued object transmission functions, particularly under decreasing SNR conditions. This approach enables us to isolate the impact of the loss function choice from the convergence behavior associated with an unknown illumination field. The reconstruction procedure is explained in more detail in Section 3 of Supplement 1, and the source code and raw data are available in [23] under open licenses. In scenarios characterized by a high SNR, no noticeable disparity in the image quality is observed between the two approaches. However, when confronted with low-SNR conditions, where the signal becomes immersed within the readout noise, the advantages of the mixed-statistics MLE loss function become clear. The reconstructed images exhibit superior quality and reveal finer details that would otherwise remain obscured without accounting for the readout noise statistics.

For a quantitative analysis and validation of our experimental findings, we generate a simulated object with phase and amplitude contrasts which we treat as the ground truth. Using numerical simulation and the ground truth object, we compute the noisy diffraction patterns assuming Poissonian photon count statistics and an additive Gaussian readout noise with a standard deviation of $\sigma = {1.5}\,\textrm{counts}$. The simulation allows for varying the illumination intensity given as the total number of photons in an otherwise fixed illumination field that approximates the experimental conditions shown in Fig. 3, thereby controlling the measurement SNR.

The achieved reconstruction quality of the object $O$ with respect to the ground truth $O_\mathrm {gt}$ can now be quantified using the correlation coefficient $C = \frac {{\left | \left \langle \overline {O}_\mathrm {gt}, O \right \rangle \right |}}{{\|O_\mathrm {gt}\| \cdot \|O\|}}$ as defined and motivated in [21] as a function of illumination intensity (Fig. 4). Here, $\overline {O}_\mathrm {gt}$ denotes the complex conjugate of $O_\mathrm {gt}$, $\left \langle \cdot, \cdot \right \rangle$ denotes the dot product, and $\|\cdot \|$ denotes the norm. The simulation confirms the same trend that we observe from the experiments: for low SNR, optimization using a mixed-statistics loss function yields significantly better reconstruction results. In the regime of high illumination intensities, all loss functions converge excellently (up to machine precision) as the readout noise becomes irrelevant, and all assumed underlying probability density functions become valid approximations. Image reconstruction using a Gaussian loss function $L_{\mathrm {Gaussian}}(\boldsymbol{\theta })=\sum _{k=1}^N\left ( X_k- I_k(\boldsymbol{\theta })\right )^2$ performs worst at low illumination intensities. However, introducing a weighting term as motivated in [24] leads to a noteworthy improvement (shown in red) using $L_{\mathrm {normMSE}}(\boldsymbol{\theta })=\sum _{k=1}^N\left ( \frac {X_k- I_k(\boldsymbol{\theta })}{\mathrm {sg}[I_k(\boldsymbol{\theta })] + \epsilon }\right )^2$, where $\epsilon =10^{-3}$ and $\mathrm {sg}[\cdot ]$ indicates a stop-gradient function. This normalized MSE loss function can be interesting in cases where $\sigma _k^2$ is impractical to obtain.

 figure: Fig. 4.

Fig. 4. Correlation between ground truth and reconstructions as a function of the total photon count in the illumination field based on numerical simulations. The dashed blue line represents reconstructions derived from intensity data where negative values have been zero-cropped, while the solid line represents reconstructions that incorporate negative values, which can arise due to the additive Gaussian noise component. (A) Amplitude and phase contrast reconstruction of the simulated object based on high-intensity diffraction patterns. For comparison, with an illumination intensity of 3.4×105 photons the reconstructions achieved with $L_\mathrm {Mixed}$ (B) and $L_\mathrm {Poisson}$ (C) are shown. The edge length of every image equals 3.5 mm.

Download Full Size | PDF

The intensity readout at a given pixel may be negative due to the additive Gaussian component in the noise statistics. This can occur in practice in an experiment via background subtraction, when areas on the camera sensor detect only low intensities. As a practical measure to keep the loss function real-valued when calculating the square root of intensity for the Poissonian loss function, negative-intensity values are customarily forced to zero [25,26]. Hence, in both simulation and experimental scenarios, we assign zero to negative-intensity values when optimizing $L_{\mathrm {Poisson}}$. However, by zero-cropping the intensity data, we may inadvertently eliminate valuable information, thereby causing a potential bias in our reconstruction results. To examine this bias and quantify the potential information contained within negative values, we also test $L_{\mathrm {Mixed}}$ with zero-cropped data, as shown in Fig. 4. The quality of reconstruction in these conditions falls between the results from optimizing $L_{\mathrm {Poisson}}$ with zero-cropped data and $L_{\mathrm {Mixed}}$ with unaltered data. This observation suggests that the enhanced reconstruction quality derived from using a mixed-statistics loss function can partially be credited to the statistical information encapsulated in the negative pixel values resulting from background subtraction.

The results presented in this study underscore the importance of considering mixed Poisson–Gaussian noise statistics in ptychographic image reconstruction. We have demonstrated, through both experimental and simulated data, that using an MLE loss function that considers this mix of noise statistics improves the reconstruction quality, particularly in low signal-to-noise ratio conditions. This enhanced performance indicates that the mixed-statistics loss function can extract more information from the measured data by accurately accounting for the underlying noise statistics. An interesting outcome of this study concerns the practice of zero-cropping negative-intensity values. We find that this practice introduces a bias into the reconstructions, highlighting the importance of preserving all statistical information in the data.

It is worth noting that some types of detectors bypass the issue of a significant Gaussian readout noise, such as high-performance photon-counting hybrid pixel detectors notably used in x ray ptychography [27,28]. In such cases, optimizing $L_{\mathrm {Poisson}}(\boldsymbol{\theta })$ can yield excellent reconstruction results. Future research could focus on studying the convergence behavior of a mixed-statistics loss function when the illumination field is jointly optimized, as we have observed that optimizing $L_{\mathrm {Mixed}}(\boldsymbol{\theta })$ occasionally leads to a less reliable convergence when dealing with a poor initial estimate for the illumination field. To ensure valid comparison and to attain image retrieval under extremely ill-posed conditions, we included the additional step of precalibrating the illumination field in this study. Such a step is typically unnecessary in well-posed ptychographic reconstructions [29] or other approaches to noise-robust phase retrieval methods [30,31].

In summary, the findings presented here could potentially propel significant advancements in the field of computational imaging, leading to improved image retrieval under challenging noise conditions. By offering a more accurate reflection of real-world ptychographic measurements, a loss function that considers mixed Poisson–Gaussian noise statistics could greatly contribute to various fields, including material science, biology, and nanotechnology, where high-quality image reconstruction under low-SNR conditions is critical. Moreover, the utility of a mixed-statistics loss function is not just limited to ptychography but extends to many computational and gradient-based imaging methods, broadening its applicability [3236].

Funding

Nederlandse Organisatie voor Wetenschappelijk Onderzoek (Perspective P16-08).

Disclosures

The authors declare no conflicts of interest.

Data availability

Raw data and source code are available in [23].

Supplemental document

See Supplement 1 for supporting content.

REFERENCES

1. J. M. Rodenburg and H. M. L. Faulkner, Appl. Phys. Lett. 85, 4795 (2004). [CrossRef]  

2. J. Rodenburg and A. Maiden, Springer Handbook of Microscopy (Springer International Publishing, 2019), p. 2.

3. K. Giewekemeyer, P. Thibault, S. Kalbfleisch, A. Beerlink, C. M. Kewish, M. Dierolf, F. Pfeiffer, and T. Salditt, Proc. Natl. Acad. Sci. U. S. A. 107, 529 (2010). [CrossRef]  

4. J. Marrison, L. Räty, P. Marriott, and P. O’Toole, Sci. Rep. 3, 2369 (2013). [CrossRef]  

5. C. C. Polo, L. Pereira, P. Mazzafera, D. N. A. Flores-Borges, J. L. S. Mayer, M. Guizar-Sicairos, M. Holler, M. Barsi-Andreeta, H. Westfahl Jr, and F. Meneau, Sci. Rep. 10, 6023 (2020). [CrossRef]  

6. D. B. Moore and J. R. Fienup, Appl. Opt. 55, 4596 (2016). [CrossRef]  

7. D. Claus, D. J. Robinson, D. G. Chetwynd, Y. Shuo, W. T. Pike, J. J. De J Toriz Garcia, and J. M. Rodenburg, J. Opt. 15, 035702 (2013). [CrossRef]  

8. P. Song, S. Jiang, H. Zhang, X. Huang, Y. Zhang, and G. Zheng, APL Photonics 4, 050802 (2019). [CrossRef]  

9. M. Du, X. Liu, A. Pelekanidis, F. Zhang, L. Loetgering, P. Konold, C. L. Porter, P. Smorenburg, K. S. E. Eikema, and S. Witte, Optica 10, 255 (2023). [CrossRef]  

10. D. Bouchet, J. Seifert, and A. P. Mosk, Opt. Lett. 46, 254 (2021). [CrossRef]  

11. H. Yang, R. N. Rutte, L. Jones, M. Simson, R. Sagawa, H. Ryll, M. Huth, T. J. Pennycook, M. L. H. Green, H. Soltau, Y. Kondo, B. G. Davis, and P. D. Nellist, Nat. Commun. 7, 12532 (2016). [CrossRef]  

12. K. Kharitonov, M. Mehrjoo, M. Ruiz-Lopez, B. Keitel, S. Kreis, M. Seyrich, M. Pop, and E. Plönjes, Opt. Express 29, 22345 (2021). [CrossRef]  

13. Z. Chen, M. Odstrcil, Y. Jiang, Y. Han, M.-H. Chiu, L.-J. Li, and D. A. Muller, Nat. Commun. 11, 2994 (2020). [CrossRef]  

14. S. M. Kay, Fundamentals of Statistical Signal Processing: Estimation Theory (Prentice-Hall, Inc., 1993).

15. P. Thibault and M. Guizar-Sicairos, New J. Phys. 14, 063004 (2012). [CrossRef]  

16. P. Godard, M. Allain, V. Chamard, and J. Rodenburg, Opt. Express 20, 25914 (2012). [CrossRef]  

17. L.-H. Yeh, J. Dong, J. Zhong, L. Tian, M. Chen, G. Tang, M. Soltanolkotabi, and L. Waller, Opt. Express 23, 33214 (2015). [CrossRef]  

18. Y. Zhang, P. Song, and Q. Dai, Opt. Express 25, 168 (2017). [CrossRef]  

19. X. Wei, H. P. Urbach, and W. M. J. Coene, Phys. Rev. A 102, 043516 (2020). [CrossRef]  

20. O. Bunk, M. Dierolf, S. Kynde, I. Johnson, O. Marti, and F. Pfeiffer, Ultramicroscopy 108, 481 (2008). [CrossRef]  

21. J. Seifert, D. Bouchet, L. Loetgering, and A. P. Mosk, OSA Continuum 4, 121 (2021). [CrossRef]  

22. M. Du, Y. S. G. Nashed, S. Kandel, D. Gürsoy, and C. Jacobsen, Sci. Adv. 6, eaay3700 (2020). [CrossRef]  

23. J. Seifert, Y. Shao, and A. P. Mosk, “Supplemental code and raw data for maximum-likelihood estimation in ptychography in the presence of poisson-gaussian noise statistics” Data publication platform of Utrecht University, January 8, (2023), https://doi.org/10.24416/UU01-4SVE6U.

24. B. Mildenhall, P. Hedman, R. Martin-Brualla, P. Srinivasan, and J. T. Barron, in 2022 IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR) (2021), pp. 16169–16178.

25. R. Cao, T. Yang, Y. Fang, C. Kuang, and X. Liu, Appl. Opt. 56, 6930 (2017). [CrossRef]  

26. S. Kandel, S. Maddali, M. Allain, S. O. Hruszkewycz, C. Jacobsen, and Y. S. G. Nashed, Opt. Express 27, 18653 (2019). [CrossRef]  

27. M. Guizar-Sicairos, I. Johnson, A. Diaz, M. Holler, P. Karvinen, H.-C. Stadler, R. Dinapoli, O. Bunk, and A. Menzel, Opt. Express 22, 14859 (2014). [CrossRef]  

28. F. Pfeiffer, Nat. Photonics 12, 9 (2018). [CrossRef]  

29. A. M. Maiden and J. M. Rodenburg, Ultramicroscopy 109, 1256 (2009). [CrossRef]  

30. C. Shen, X. Bao, J. Tan, S. Liu, and Z. Liu, Opt. Express 25, 16235 (2017). [CrossRef]  

31. X. Chang, C. Shen, S. Liu, D. Zheng, S. Wang, C. Yang, N. E. Huang, and L. Bian, Opt. Lett. 48, 4161 (2023). [CrossRef]  

32. N. Antipa, G. Kuo, R. Heckel, B. Mildenhall, E. Bostan, R. Ng, and L. Waller, Optica 5, 1 (2018). [CrossRef]  

33. M. O’Toole, D. B. Lindell, and G. Wetzstein, Nature 555, 338 (2018). [CrossRef]  

34. R. Raskar, A. Agrawal, C. A. Wilson, and A. Veeraraghavan, ACM Trans. Graph. 27, 1 (2008). [CrossRef]  

35. S.-H. Hong, J.-S. Jang, and B. Javidi, Opt. Express 12, 483 (2004). [CrossRef]  

36. A. Liutkus, D. Martina, S. Popoff, G. Chardon, O. Katz, G. Lerosey, S. Gigan, L. Daudet, and I. Carron, Sci. Rep. 4, 5552 (2014). [CrossRef]  

Supplementary Material (1)

NameDescription
Supplement 1       Supplemental document

Data availability

Raw data and source code are available in [23].

23. J. Seifert, Y. Shao, and A. P. Mosk, “Supplemental code and raw data for maximum-likelihood estimation in ptychography in the presence of poisson-gaussian noise statistics” Data publication platform of Utrecht University, January 8, (2023), https://doi.org/10.24416/UU01-4SVE6U.

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (4)

Fig. 1.
Fig. 1. Comparative visualization of a grayscale image distorted by different noise types. Panel (A) depicts an image with simulated Poissonian noise, while panel (B) illustrates the effects of simulated mixed Poisson–Gaussian noise resulting from additive readout noise. Inlaid values indicate peak signal-to-noise ratio (PSNR) with respect to the ground truth.
Fig. 2.
Fig. 2. Schematic drawing of our ptychography setup used in the experiment and for numerical simulations. A 500-µm pinhole is illuminated and relayed onto the object using a 2-lens system. The object is shifted laterally through the beam using an XY stage. A CMOS camera records the diffraction intensities at a distance of 38 mm downstream of the object.
Fig. 3.
Fig. 3. Comparison of image reconstruction qualities from ptychographic datasets with decreasing signal-to-noise ratio (SNR) from left to right. Top row: Ptychographic datasets at different camera exposures with total illumination intensity, respectively. Second row: Reconstructions using a loss function assuming solely Poissonian noise statistics. Third row: Reconstructions using the loss function defined by Eq. (2), which incorporates mixed Poisson–Gaussian noise statistics. The edge length of every image equals 3.5 mm.
Fig. 4.
Fig. 4. Correlation between ground truth and reconstructions as a function of the total photon count in the illumination field based on numerical simulations. The dashed blue line represents reconstructions derived from intensity data where negative values have been zero-cropped, while the solid line represents reconstructions that incorporate negative values, which can arise due to the additive Gaussian noise component. (A) Amplitude and phase contrast reconstruction of the simulated object based on high-intensity diffraction patterns. For comparison, with an illumination intensity of 3.4×105 photons the reconstructions achieved with $L_\mathrm {Mixed}$ (B) and $L_\mathrm {Poisson}$ (C) are shown. The edge length of every image equals 3.5 mm.

Equations (2)

Equations on this page are rendered with MathJax. Learn more.

L P o i s s o n ( θ ) = k = 1 N ( X k I k ( θ ) ) 2 ,
L M i x e d ( θ ) = k = 1 N ( ln [ I k ( θ ) + σ k 2 ] + [ X k I k ( θ ) ] 2 I k ( θ ) + σ k 2 ) .
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.