Deep Ptych: Subsampled Fourier Ptychography using Generative Priors

TL;DR

Deep Ptych introduces a generative prior-based method for Fourier ptychography, significantly enhancing reconstruction quality and noise robustness with fewer samples, and allowing exploration beyond the generative model's range.

Contribution

The paper presents a novel deep learning framework that regularizes Fourier ptychography using generative models, improving reconstruction quality and robustness.

Findings

Outperforms existing techniques in reconstruction quality

Requires fewer samples for accurate results

Enhances robustness against noise

Abstract

This paper proposes a novel framework to regularize the highly ill-posed and non-linear Fourier ptychography problem using generative models. We demonstrate experimentally that our proposed algorithm, Deep Ptych, outperforms the existing Fourier ptychography techniques, in terms of quality of reconstruction and robustness against noise, using far fewer samples. We further modify the proposed approach to allow the generative model to explore solutions outside the range, leading to improved performance.

Figures

Tables

Table 1. Table 1: PSNR (dB) results for MNIST and CelebA.

	MNIST				CelebA
	Subsampling		Noise		Subsampling		Noise
	1%	3%	1%	10%	1%	3%	1%	10%
IERA	10.59	11.72	14.15	9.59	5.89	7.29	12.36	11.25
CoPRAM	14.17	19.16	16.63	6.45	13.95	23.85	17.53	9.34
DP	24.70	24.80	25.31	19.27	22.54	22.79	23.57	21.21
DP+	29.78	39.39	34.83	15.20	27.99	38.72	32.08	14.00
DP+(TV)	-	-	35.98	16.82	-	-	33.61	16.46

Equations

$${\mathbf{y}}_{\ell }=|{\mathcal{A}}_{\ell }(\mathbf{x})|+{\mathbf{n}}_{\ell },\text{for}\ell =1,2,...,L,$$

$${\mathbf{y}}_{\ell }=\sum _{\ell =1}^L|{\mathcal{A}}_{\ell }(\mathbf{x})|+{\mathbf{n}}_{\ell },$$

$$\text{Subsampling Ratio (\%)}=\frac{\text{Fraction of samples retained (f)}\times 100}{\text{Total observed samples (nL)}}.$$

$$\hat{\mathbf{x}}:=\underset{\begin{array}{c}\mathbf{x}\in \text{Range}(G)\end{array}}{\mathrm{argmin}}\sum _{\ell =1}^L\Vert {\mathbf{y}}_{\ell }-|{\mathcal{A}}_{\ell }(x)|{\Vert }_2^2,$$

$$\hat{\mathbf{z}}=\underset{\mathbf{z}\in {\mathbb{R}}^k}{\mathrm{argmin}}\sum _{\ell =1}^L\Vert {\mathbf{y}}_{\ell }-|{\mathcal{A}}_{\ell }(G(\mathbf{z}))|{\Vert }_2^2.$$

$$(\hat{\mathbf{x}},\hat{\mathbf{z}})=\underset{\mathbf{z}\in {\mathbb{R}}^k,\mathbf{x}\in {\mathbb{R}}^n}{\mathrm{argmin}}\sum _{\ell =1}^L\Vert {\mathbf{y}}_{\ell }-|{\mathcal{A}}_{\ell }(\mathbf{x})|{\Vert }_2^2+\lambda \Vert \mathbf{x}-G(\mathbf{z}){\Vert }_2^2.$$

Full text

1 Introduction

Resolution loss in long distance imaging can primarily be attributed to the diffraction blur, that is caused by limited aperture of the imaging system [1]. To mitigate the effects of the diffraction blur, recently an emerging computational imaging technique known as Fourier Ptychography (FP) has shown promising results [2, 3]. The objective of FP is to recover a high-resolution image from multiple diffraction-limited low-resolution images. In this paper, we consider recovering the signal ${\boldsymbol{x}}\in\mathbb{C}^{n}$ captured via forward acquisition model of FP, given by:

[TABLE]

where ${\boldsymbol{y}}_{\ell}\in\mathbb{R}^{n}$ is diffraction-limited image corresponding to $\ell^{th}$ camera, ${\mathcal{A}}_{\ell}:\mathbb{C}^{n}\rightarrow\mathbb{C}^{n}$ is the linear operator representing the forward acquisition model, and ${\boldsymbol{n}}_{\ell}\in\mathbb{R}^{n}$ denotes noise perturbation. For $\ell^{th}$ camera, the linear operator ${\mathcal{A}}_{\ell}$ has the form ${\mathcal{F}}^{-1}{\mathcal{P}}_{\ell}{\circ}{\mathcal{F}}$, where ${\mathcal{F}}$ denotes 2D Fourier transform, ${\mathcal{P}}_{\ell}$ is a pupil mask that acts as a bandpass filter in the Fourier domain, and $\circ$ represents the Hadamard product (details are provided in subsequent sections). Specifically, FP works by iteratively stitching together a sequence of frequency limited low-resolution images ${\boldsymbol{y}}_{\ell}$ in Fourier domain to recover the high-resolution true image ${\boldsymbol{x}}$. Since optical sensors can measure only the magnitude of the signal [4], phase information is lost during the acquisition process — making the FP problem highly ill-posed.

Traditional approaches to overcome the ill-posedness of FP are iterative phase retrieval algorithms [5]. Generally, phase retrieval algorithms fall into two categories. The first approach is to introduce redundancy into the measurement system, where more measurements are taken than the dimension of true signal usually in the form of oversampled Fourier transform [6], short-time Fourier transform [7], structured illuminations [8], etc. The second approach is to exploit the structural assumptions on the true signal as a prior, such as sparsity [9] or non-negativity [10]. These prior based phase retrieval approaches have recently attained much significance for the purpose of reducing number of measurements as the acquisition of measurements is usually expensive and time-consuming, especially at high resolutions.

FP recovery algorithms that exploit redundant information perform well provided the set of low-resolution images ${\boldsymbol{y}}_{\ell}$ have high overlapping frequency bands (typically $60\%$ overlap) in the Fourier domain [11]. The reconstruction quality of these algorithms quickly degrades with the decrease in overlap. This high overlap requirement for faithful reconstruction of the true image ${\boldsymbol{x}}$ results in long acquisition times and high computational cost [2].

To reduce the computational cost, algorithms that exploit prior information about true signal have recently been explored for FP. Among them, priors learned from large datasets by utilizing the power of deep neural networks have shown promising results for low overlap case, thus reducing computational cost [12, 13, 14]. Specifically, these deep learning based approaches invert the forward acquisition model of FP via end-to-end training of deep neural networks in a supervised manner [15]. However, even a slight change in the noise level ${\boldsymbol{n}}_{\ell}$ or parameters of FP measurement model, such as aperture diameter or overlap would require costly retraining of these models. Apart from reducing overlap, an alternative approach to reduce computational cost for FP, is to design realistic subsampling strategies by exploiting some prior information about the true signal. This alternative has not yet been investigated by the deep learning based FP approaches. Recently, Jagatap et al. [16] explored this alternative by leveraging sparsity prior. The resulting algorithm (CoPRAM - Compressive Phase Retrieval using Alternative Minimization) significantly reduces the number of samples required for faithful reconstruction of the true signal ${\boldsymbol{x}}$ in FP setup. However, it has been shown in [17] that sparsity priors often fail to capture the complex structure that many natural signals exhibit resulting in unrealistic signals also fitting the sparse prior modelling assumption.

Recently, neural network based implicit generative models such as Generative Adversarial Networks (GANs) [18] and Variational Autoencoders (VAEs) [19] have been quite successful in modelling complex data distributions especially that of images. In these models, the generative part $(G)$, learns a mapping from low dimensional latent space ${\boldsymbol{z}}\in\mathbb{R}^{k}$ to a high dimensional sample space $G({\boldsymbol{z}})\in\mathbb{R}^{n}$ where $k\ll n$. During training, these generative models are encouraged to produce samples that resemble with that of training data $\mathcal{X}$. A well-trained generator, given by deterministic function $G:\mathbb{R}^{k}\rightarrow\mathbb{R}^{n}$ with a distribution $P_{Z}$ over ${\boldsymbol{z}}$ (usually random normal or uniform), is therefore capable of generating fake data indistinguishable from the real data it has been trained on. Due to their power to effectively model natural image distributions, generative priors have recently been introduced to solve ill-posed inverse problems like compressed sensing [20], phase retrieval [21, 22], blind deconvolution [23], etc. Notably, these generative prior based approaches, have been shown to improve over sparsity-based approaches, thus advancing the state of the art in several fields [24].

1.1 Our Contributions

In this paper, we leverage the expressive power of generative models to regularize FP problem. To the best of our knowledge, this is the first work that explores the generative models for FP setup. In particular, we make the following contributions:

We combine the expressive power of deep generative models (GANs and VAEs) with non-linear inverse problem of FP for the first time. We show that the resultant problem can be solved effectively using conventional gradient descent algorithm yielding promising results. Moreover, unlike deep learning based approaches [12, 14], our approach does not require expensive retraining in case of slight changes in parameters such as subsampling ratio or noise level. 2. 2.

We show experimentally that the proposed approach produces impressive results for low subsampling ratios and high noise. 3. 3.

Finally, to further improve the performance of the proposed scheme, especially for the rich image datasets not well learned by the generator, we present a modification in the loss function that allows more flexible reconstructions.

2 Formulation

In this section, we describe the forward acquisition model of FP for long-distance imaging following the setup of [2], as shown in Fig. 1. In long-distance imaging, resolution of the image is limited by diffraction of light at the limited aperture of the imaging system. To emulate capturing the scene with a large synthetic aperture, multiple cameras ($L$) are deployed, usually in a square grid. Such an arrangement is called a coherent camera array. In long-distance imaging, coherent camera array is usually placed in far field of the object (${\boldsymbol{x}}$), and the Fraunhofer approximation is satisfied. Under this approximation, illumination field emerging from the object when intercepted by the camera array can be treated as 2D Fourier transform (${\mathcal{F}}$) of the image at the object plane [25]. Each camera in the array have a limited aperture denoted by ${\mathcal{P}}_{\ell}$ ($\ell=1,2,...,L$) that acts as a bandpass filter covering different parts of Fourier domain image as shown in Fig. 1. The bandlimited signal is then focused to an image plane where this propagation to the image plane can be modelled as an inverse Fourier transform (${\mathcal{F}}^{-1}$). After that, a subsampling operator (${\mathcal{M}}_{\ell}(\cdot)$), is applied to effectively reduce the number of measurements. Subsequently, the complex spatial domain image is captured by an optical sensor that measures only the magnitude while discarding the phase information. Mathematically, for a total of $L$ cameras, observations ${\boldsymbol{y}}_{\ell}$ can be modelled as

[TABLE]

where ${\mathcal{A}}_{\ell}={\mathcal{M}}_{\ell}{\mathcal{F}}^{-1}{\mathcal{P}}_{\ell}\circ{\mathcal{F}}$ is the measurement model prior to optical sensor acquisition step, $\circ$ denotes the Hadamard product and ${\mathcal{M}}_{\ell}$ is the subsampling operator. Subsampling operator when applied to an image, randomly picks a fraction of samples $(f)$ discarding the others [16]. We define the subsampling ratio as the fraction of samples retained by ${\mathcal{M}}_{\ell}$ divided by the total number of observed samples i.e.

[TABLE]

3 Deep Ptych

In this section, we formally introduce our proposed approach Deep Ptych. We assume that the true image $\boldsymbol{x}\in{\mathbb{R}}^{n}$ is a member of some structured class of images, denoted by $\mathcal{X}$. For example, $\mathcal{X}$ may be set of face images or digits. A generative model is trained on a representative sample set from class $\mathcal{X}$. Given low-dimensional input vector $\boldsymbol{z}\in\mathbb{R}^{k}$, the generator $G$ after training, generates new samples $G(\boldsymbol{z})$ similar to the representative samples of class $\mathcal{X}$. Once trained, the weights of the generator are kept fixed (pretrained generator). To recover the estimate of true image $\hat{{\boldsymbol{x}}}$ from FP measurements $\boldsymbol{y}_{\ell}$ for $\ell=1,2,...,L,$ in (2), we propose minimizing the following objective function

[TABLE]

where Range($G$) is the set of all the images that can be generated by pretrained $G$. In other words, we want to find an image $\hat{{\boldsymbol{x}}}$ that is closest to the observed measurements and lies within the range of the generator. Ideally, the range of the generator must comprise of only the samples drawn from the distribution of the images it has been trained on, i.e., $\exists\,{\boldsymbol{z}}$ such that ${\boldsymbol{x}}=G(\boldsymbol{z})$. It must be noted that constraining the solution ${\boldsymbol{\hat{x}}}$ to lie inside the generator range will force the solution to be a member of image class $\mathcal{X}$ thus, implicitly reducing the solution ambiguities associated with FP [5]. The minimization program in (3) can be equivalently formulated in the lower dimensional, latent representation ${\boldsymbol{z}}$ as follows:

[TABLE]

This optimization program modifies the latent representation vector $\boldsymbol{z}$ such that the generator generates an image $G({\boldsymbol{z}})$ that is consistent with (4). Since the latent representation ${\boldsymbol{z}}$ is much smaller in dimension as compared to the ambient dimension of the true image ${\boldsymbol{x}}$, it not only reduces the number of unknowns in the FP problem but also allows an efficient implementation of gradients using backpropagation through the generators. Since optimization program in (4) is non-convex and non-linear owing to the modulus operator and non-linear generative model, we resort to the gradient descent algorithm to find an optimal $\hat{{\boldsymbol{z}}}$. The estimated image $\hat{{\boldsymbol{x}}}$ is acquired by a forward pass of the $\hat{{\boldsymbol{z}}}$ through the generator $G$. Mathematically, $\hat{{\boldsymbol{x}}}=G(\hat{{\boldsymbol{z}}})$. The gradient descent scheme, namely, Deep Ptych is formally given in Algorithm 1, where the objective function in (4) is denoted by $\mathcal{L}({\boldsymbol{z}})$, that is, $\mathcal{L}({\boldsymbol{z}}):=\sum_{\ell=1}^{L}\|{\boldsymbol{y}}_{\ell}-|\mathcal{A_{\ell}}({G}({\boldsymbol{z}}))|\|^{2}_{2}$, and $\nabla\mathcal{L}({\boldsymbol{z}})$ denotes its gradient w.r.t. ${\boldsymbol{z}}$.

3.1 Beyond the range of generator: Deep Ptych+

The shortcoming of the generative network in effectively learning the distribution of the image dataset, especially more complex ones, hampers its ability to reliably reproduce a given new (test) sample from its range [26]. To address this shortcoming, we allow the reconstructed image to deviate a bit from the range of the generator if it improves the measurement loss. To achieve this, we propose solving the following modified version of the optimization program in (4)

[TABLE]

The first term in the objective favors an image ${\boldsymbol{x}}$ with smaller measurement loss, while the second term ensures that in doing so such an ${\boldsymbol{x}}$ should not deviate too far away from the range of the pretrained generative network. The free parameter $\lambda\geq 0$ controls the degree to which we are willing to deviate away from the range of the generator, for example, one might choose a smaller $\lambda$ when the generator is known to not reliably reproduce the image samples from the dataset.

To minimize (5), we use an alternating gradient descent scheme that takes gradient steps in one of the unknowns ${\boldsymbol{x}}$ or ${\boldsymbol{z}}$ while keeping the other fixed. This enables us to land at the local minima $(\hat{{\boldsymbol{x}}},\hat{{\boldsymbol{z}}})$ of the objective in (5). We take $\hat{{\boldsymbol{x}}}$ to be the reconstructed image, which generally improves over the reconstructions from Deep Pytch algorithm above. We refer to this modified approach as Deep Ptych+.

4 Experimental results

We now provide simulation experiments to evaluate the performance of Deep Ptych and Deep Ptych+. To quantitatively evaluate the performance of our algorithm, we use two metrics Peak Signal to Noise Ratio (PSNR) and Structural Similarity Index Measure (SSIM). In all our experiments, we report results on a held out test set, unseen by the generative model during training.

4.1 Datasets and Generative models

We use three datasets (MNIST [27], CelebA (64$\times$64) [28], and high-resolution CelebA (128$\times$128) [29]) and two generative models (VAE and GAN). This provides some evidence that our approach is generalizable to many types of datasets and generative models.

MNIST with VAE: MNIST consists of 28$\times$28 grayscale images of handwritten digits, with 60,000 training and 10,000 test examples. We super-resolve original MNIST dataset to size 56$\times$56 using super-resolution method proposed in [30], that gives promising results without requiring external training data. We train VAE, having same architecture as proposed in [23] for image datasets, on this super-resolved dataset with size of latent dimension set to 50, batch size of 128, $\lambda$ set to $2\times 10^{-4}$, and learning rate of $5\times 10^{-4}$ with Adam optimizer.

CelebA with DCGAN: We use aligned and cropped version of CelebA dataset having more than 200,000 RGB face images of size 64$\times$64$\times$3. Each pixel value is scaled between $[-1,1]$. Deep convolutional generative adversarial network (DCGAN) is trained on this rescaled dataset with architecture as proposed in [31] having latent dimension size set to 100, batch size of 64, $\beta_{1}$ set to 0.5, $\lambda$ set to $1.4\times 10^{-4}$, and learning rate of $2\times 10^{-4}$ using Adam optimizer.

CelebA (128$\times$128) with Progressive GAN [29]: To show the dependence of Deep Ptych on the expressive power of generator, we use a powerful generative model named progressive GAN [29], trained on 128$\times$128 face images. We use the originally proposed architecture having the latent dimension set to 512.

4.2 Baseline methods

We consider IERA (Iterative Error Reduction Algorithm) [2] and CoPRAM [16] as baseline methods. IERA solves FP problem by alternatively enforcing spatial and Fourier domain constraints whereas CoPRAM exploits the underlying sparse structure of true signal for faithful reconstruction at low subsampling ratios. For CoPRAM, we assume sparsity for MNIST and CelebA in spatial and Fourier domains respectively. We use the default algorithmic parameters for IERA and CoPRAM, unless stated otherwise.

4.3 Subsampling results

We observed that CoPRAM does not work well with additive noise therefore, to carry out a fair comparison, we perform all subsampling experiments without noise. For MNIST, we run gradient descent algorithm using Adam optimizer with learning rate of 0.05 for 2000 steps. We set aperture diameter to 15 pixels and overlap to 65$\%$ with camera array consisting of $81$ cameras (square grid of $9\times 9$). We use same parameters for CelebA, except learning rate and aperture diameter, that are set to 0.01 and 16 pixels respectively.

We observe that Deep Ptych and Deep Ptych+ are able to attain much higher SSIM and PSNR values at low subsampling ratios as shown in Fig. 2 and Table 1 respectively. Qualitative results at $2\%$ subsampling ratio show faithful reconstruction of Deep Ptych and Deep Ptych+ as compared to IERA and CoPRAM results. Since output of Deep Ptych is constrained to lie in generator range, increasing subsampling ratio results in its performance saturation. CoPRAM has no such limitation so its performance eventually surpasses Deep ptych at higher subsampling ratios. We show in Fig. 4 that for a more expressive generator this range constraint becomes less effective, resulting in visually appealing reconstructions for Deep Ptych. Deep Ptych+, on the other hand, outperforms IERA and CoPRAM at all subsampling ratios.

4.4 Noise results

For all experiments in this section a noise of $1\%$, for image scaled between 0 to 1, translates to Gaussian noise with zero mean and a standard deviation of $0.01$. All parameter values are set to same as in subsampling experiments with subsampling ratio set to $10\%$. Quantitative results in Fig. 2 and Table 1 manifest that Deep Ptych exhibits far superior performance than baseline approaches for high noise levels. On the other hand, performance of Deep Ptych+ is superior to IERA and CoPRAM, but degrades quickly by increasing noise percentage. This is because its output is not constrained to lie in generator range (that spans the set of representative samples of noiseless training data). We find that by adding total variation (TV) regularization (with weighting factor of $10^{-4}$) to Deep Ptych+ loss function, the performance is improved (generally $0.5-2$dB gain in PSNR) at all noise levels. We dubbed this modified Deep Ptych+ as Deep Ptych+ (TV). Qualitatively, as depicted in Fig. 3, the proposed approaches outperform baseline methods in terms of quality of reconstruction at $2.5\%$ noise.

5 Conclusion

To conclude, we demonstrated the effectiveness of integrating deep generative priors with Fourier ptychography problem. We showed that the proposed approach is effective at low subsampling ratios and is highly robust to noise. We further refined the proposed algorithm to allow the generative model to explore solutions outside the range, leading to improved performance.

Acknowledgement

We gratefully acknowledge the support of the NVIDIA Corporation for the donation of NVIDIA TITAN Xp GPU for this research.