Smart quantum statistical imaging beyond the Abbe-Rayleigh criterion

The wave nature of light imposes limits on the resolution of optical imaging systems. For over a century, the Abbe-Rayleigh criterion has been utilized to assess the spatial resolution limits of imaging instruments. Recently, there has been interest in using spatial projective measurements to enhance the resolution of imaging systems. Unfortunately, these schemes require a priori information regarding the coherence properties of “unknown” light beams and impose stringent alignment conditions. Here, we introduce a smart quantum camera for superresolving imaging that exploits the self-learning features of artificial intelligence to identify the statistical fluctuations of unknown mixtures of light sources at each pixel. This is achieved through a universal quantum model that enables the design of artificial neural networks for the identification of photon fluctuations. Our protocol overcomes limitations of existing superresolution schemes based on spatial mode projections, and consequently provides alternative methods for microscopy, remote sensing, and astronomy.

The spatial resolution of optical imaging systems is established by the diffraction of photons and the noise associated with their quantum fluctuations [1][2][3][4][5].For over a century, the Abbe-Rayleigh criterion has been used to assess the diffraction-limited resolution of optical instruments [3,6].At a more fundamental level, the ultimate resolution of optical instruments is established by the laws of quantum physics through the Heisenberg uncertainty principle [7][8][9].In classical optics, the Abbe-Rayleigh resolution criterion stipulates that an imaging system cannot resolve spatial features smaller than λ/2NA.In this case, λ represents the wavelength of the illumination field, and NA describes numerical aperture of the optical instrument [1][2][3]10].Given the implications that overcoming the Abbe-Rayleigh resolution limit has for multiple applications, such as, microscopy, remote sensing, and astronomy [3,[10][11][12], there has been an enormous interest in improving the spatial resolution of optical systems [13][14][15].So far, optical superresolution has been achieved through spatial decomposition of eigenmodes [14,16,17].These conventional schemes rely on spatial projective measurements to pick up phase information that is used to boost spatial resolution of optical instruments [14,[18][19][20][21][22].
The quantum statistical fluctuations of photons establish the nature of light sources [30][31][32][33][34].As such, these fundamental properties are not affected by the spatial resolution of an optical instrument [34].Here, we demonstrate that measurements of the quantum statistical properties of a light field enable imaging beyond the Abbe-Rayleigh resolution limit.This is performed by exploiting the self-learning features of artificial intelligence to identify the statistical fluctuations of photon mixtures [30].More specifically, we demonstrate a smart quantum camera with the capability to identify photon statistics at each pixel.For this purpose, we introduce a universal quantum model that describes the photon statistics produced by the scattering of an arbitrary number of light sources.This model is used to design and train artificial neural networks for the identification of light sources.

a b
FIG. 1. Conceptual illustration and schematic of our experimental setup to demonstrate superresolving imaging.The illustration in a depicts a scenario where diffraction limits the resolution of an optical instrument for remote imaging.In our protocol, an artificial neural network enables the identification of the photon statistics that characterize the point sources that constitute a target object.In this case, the point sources emit either coherent or thermal photons.Remarkably, the neural network is capable of identifying the corresponding photon fluctuations and their combinations, for example coherent-thermal (CT1, CT2), thermal-thermal (TT) and coherent-thermal-thermal (CTT).This capability allows us to boost the spatial resolution of optical instruments beyond the Abbe-Rayleigh resolution limit.The experimental setup in b is designed to generate two independent thermal and one coherent light sources.The three sources are produced from a continuous-wave (CW) laser at 633 nm.The CW laser beam is divided by two beam splitters (BS) to generate three spatial modes, two of which are then passed through rotating ground glass (RGG) disks to produce two independent thermal light beams.The three light sources, with different photon statistics, are attenuated using neutral density (ND) filters and then combined to mimic a remote object such as the one shown in the inset of b.This setup enables us to generate multiple sources with tunable statistical properties.
The generated target beam is then imaged onto a digital micro-mirror device (DMD) that we use to perform raster scanning.The photons reflected off the DMD are collected and measured by a single-photon detector.Our protocol is formalized by performing photon-number-resolving detection [30].The characteristic quantum fluctuations of each light source are identified by an artificial neural network.This information is then used to produce a high-resolution image of the object beyond the diffraction limit.
The conceptual schematic behind our experiment is depicted in Fig. 1a.This camera utilizes an artificial neural network to identify the photon statistics of each point source that constitutes a target object.The description of the photon statistics produced by the scattering of an arbitrary number of light sources is achieved through a general model that relies on the quantum theory of optical coherence introduced by Sudarshan and Glauber [34][35][36].We use this model to design and train a neural network capable of identifying light sources at each pixel of our camera.This unique feature is achieved by performing photon-number-resolving detection [30].The sensitivity of this camera is limited by the photon fluctuations, as stipulated by the Heisenberg uncertainty principle, and not by the Abbe-Rayleigh resolution limit [5,34].
In general, realistic imaging instruments deal with the detection of multiple light sources.These sources can be either distinguishable or indistinguishable [3,34].The combination of indistinguishable sources can be repre-sented by either coherent or incoherent superpositions of light sources characterized by Poissonian (coherent) or super-Poissonionan (thermal) statistics [34].In our model, we first consider the indistinguishable detection of N coherent and M thermal sources.For this purpose, we make use of the P-function P coh (γ) = δ 2 (γ −α k ) to model the contributions from the kth coherent source with the corresponding complex amplitude α k [35,36].The total complex amplitude associated to the superposition of an arbitrary number of light sources is given by α tot = N k=1 α k .In addition, the P-function for the lth thermal source, with the corresponding mean photon numbers ml , is defined as P th (γ) = (π ml ) −1 exp (−|γ| 2 / ml ).The total number of photons attributed to the M number of thermal sources is defined as m tot = M l=1 ml .These quantities allow us to calculate the P-function for the multisource system as (1) ω ij FIG. 2. The scheme of the two-layer neural network used to identify the photon statistics produced by a combination of three sources is shown in a.The computational model consists of an input layer, a hidden layer of sigmoid neurons, and a Softmax output layer.The training of our neural network through Eqs. 13 and 3 enables the efficient identification of five classes of photon statistics.Each class is characterized by a g (2) function, which is defined by a specific combination of light sources [30].In our experiment, these classes correspond to the characteristic photon statistics produced by coherent or thermal light sources and their combinations.For example, coherent-thermal, thermal-thermal, or coherent-thermal-thermal.The figure in b shows the performance of our neural network as a function of the number of data samples used each time in the testing process.The classification accuracy for the five possible complex classes of light is 80% with 100 data points.Remarkably, the performance of the neural network increases to approximately 95% when we use 3500 data points in each test sample.
This approach enables the analytical description of the photon-number distribution p th-coh (n) associated to the detection of an arbitrary number of indistinguishable light sources.This is calculated as p th-coh (n) = n| ρth-coh |n , where ρ th-coh = P th-coh (γ) |γ γ|d 2 γ.After algebraic manipulation (see Supplementary Information), we obtain the following photon-number distribution where Γ(z) and 1 F 1 (a; b; z) are the Euler gamma and the Kummer confluent hypergeometric functions, respectively.This probability function enables the general description of the photon statistics produced by any indistinguishable combination of light sources.Thus, the photon distribution produced by the distinguishable detection of N light sources can be simply obtained by performing a discrete convolution of Eq. 13 as (3) The combination of Eq. 13 and Eq. 3 allows the classification of photon-number distributions for any combination of light sources.
We demonstrate our proof-of-principle quantum camera using the experimental setup shown in Fig. 1b.For this purpose, we use a continuous-wave laser at 633nm to produce either coherent, or incoherent superpositions of distinguishable, indistinguishable, or partially distinguishable light sources.In this case, the combination of photon sources, with tunable statistical fluctuations, acts as our target object.Then, we image our target object onto a digital micro-mirror device (DMD) that is used to implement raster scanning.This is implemented by selectively turning on and off groups of pixels in our DMD.The light reflected off the DMD is measured by a single-photon detector that allows us to perform photon-number-resolving detection.This is implemented through the technique described in ref. [30].
The equations above allow us to implement a multilayer feed-forward network for the identification of the quantum photon fluctuations of the point sources of a target object.The structure of the network consists of a group of interconnected neurons arranged in layers.Here, the information flows only in one direction, from input to output [37,38].As indicated in Fig. 2a, our network comprises two layers, with ten sigmoid neurons in the hidden layer (green neurons) and five softmax neurons in the output layer (orange neurons).In this case, the input features represent the probabilities of detecting n photons at a specific pixel, p(n), whereas the neurons in the last layer correspond to the classes to be identified.The input vector is then defined by twenty-one features corresponding to n=0,1,...,20.In our experiment, we define five classes that we label as: coherent-thermal (CT), thermal-thermal (TT), coherent-thermal-thermal (CTT), coherent (C), and thermal (T).If the brightness of the experiment remains constant, these classes can be directly defined through the photon-number distribution described by Eqs. 13 and 3.However, if the brightness of the sources is modified, the classes can be defined through the g (2) = 1 + (∆n) 2 − n / n 2 , which is intensity-independent [30,31].The parameters in the g (2) function can also be calculated from Eqs. 13 and 3.It is important to mention that the output neurons provide a probability distribution over the predicted classes [39,40].The training details of our neural networks can be found in the Methods section.
We test the performance of our neural network through the classification of a complex mixture of photons produced by the combination of one coherent with two thermal light sources.The accuracy of our trained neural network is reported in Fig. 2b.In our setup, the three partially overlapping sources form five classes of light with different mean photon numbers and photon statistics.We exploit the functionality of our artificial neural network to identify the underlying quantum fluctuations that characterize each kind of light.We calculate the accuracy as the ratio of true positive and true negative to the total of input samples during the testing phase.Fig. 2b shows the overall accuracy as a function of the number of data points used to build the probability distributions for the identification of the multiple light sources using a supervised neural network.The classification accuracy for the mixture of three light sources is 80% with 100 photon-number-resolving measurements.The performance of the neural networks increases to approximately 95% when we use 3500 data points to generate probability distributions.
The performance of our protocol for light identification can be understood through the distribution of light sources in the probability space shown in Fig. 3.Here we show the projection of the feature space on the plane defined by the probabilities p(0), p(1), and p(2) for different number of data points.Each point is obtained from an experimental probability distribution.As illustrated in Fig. 3a, the distributions associated to the multiple sources obtained for 10 data points are confined to a small region of the feature space.This condition makes extremely hard the identification of light sources with 10 sets of measurements.A similar situation can be observed for the distribution in Fig. 3b that was generated using 100 data points.As shown in panel Fig. 3c, the separations in the distributions produced with 1000 data points occupy different regions, although brown and black points keep closely intertwined.These conditions enable one to identify multiple light sources.Finally, the separated distributions obtained with 10000 data points in Fig. 3d enable efficient identification of light sources.These probability space diagrams explain the performances reported in Fig. 2.An interesting feature of Fig. 3 is the fact that the distributions in the probability space are linearly separable.As demonstrated in Fig. 4, the identification of the quantum photon fluctuations at each pixel of our camera enables us to demonstrate superresolving imaging.In our experiment we prepared each source to have a mean photon number between 1 and 1.5 for the brightest pixel.The raster-scan image of a target object composed of multiple partially distinguishable sources in Fig. 4a illustrates the performance of conventional imaging protocols limited by diffraction [4,[6][7][8].In this case, it is practically impossible to identify the multiple sources that constitute the target object.Remarkably, as shown in Fig. 4b, our protocol provides a dramatic improvement of the spatial resolution of the imaging system.In this case, it becomes clear the presence of the three emitters that form the remote object.The estimation of separations among light sources is estimated through a fit over the classified pixel-by-pixel image.Additional details can be found in the Methods section.In Figs.4c and d, we demonstrate the robustness of our protocol by performing superresolving imaging for a different configuration of light sources.In this case, two small sources are located inside the point-spread function of a third light source.As shown in Fig. 4c, the Abbe-Rayleigh limit forbids the identification of light sources.However, we demonstrate substantial improvement of spatial resolution in Fig. 4d.The plots in Figs.4e and f correspond to the inferred spatial distributions based on the experimental pixel-by-pixel imaging used to produce Figs.4b  and d.The insets in Figs.4e and f show photon-number probability distributions for three pixels.The theoretical photon-number distributions in Fig. 4e and f are obtained through a procedure of least square regression [41].
Here the least squares difference between the measured and theoretical probability distribution was minimized for 0 ≤ n ≤ 6.The sources were assumed to be partially distinguishable allowing the theoretical distribution to be defined by Eqs. 13 and Eq. 3. The combined mean photon numbers of each source generated for the fit totals the measured mean photon number (see Methods section).Our scheme enables the use of the photon-number distributions or their corresponding g (2) to characterize light sources.This allows us to determine each pixel's corresponding statistics, regardless of the mean photon numbers of the sources in the detected field [30,31].FIG. 5. Comparison between the spatial resolution of our camera and direct imaging.Here the distance is normalized by the beam radius for easy identification of the Abbe-Rayleigh limit.As shown in a, the red line is the result of a Monte-Carlo simulation for traditional intensity based direct imaging.The plateau is the area where the algorithm becomes unstable.The dotted blue line represents the limit for our supperresolving imaging method, where perfect classification of each pixel is assumed.The blue dots represent the experimental data collected with our camera for superresolving imaging.The experimental points demonstrate the potential of our technique for identifying spatial features beyond the Abbe-Rayleigh resolution criterion.The first row in the panels from b to d shows the reconstructed spatial profiles obtained through direct imaging whereas the second row shows the superresolving images obtained with our technique.The panel in b shows the spatial profiles for the experimental point i).This corresponds to the experimental detection of two sources with the largest separation.The spatial profiles in c correspond to the experimental point labeled as ii).Finally, the panel in d shows the spatial distributions for the experimental point with the smallest separation, this is labeled as iii).
We now provide a quantitative characterization of our superresolving imaging scheme based on the identification of photon statistics.We demonstrate that our smart camera for superresolving imaging can capture small spatial features that surpass the resolution capabilities of conventional schemes for direct imaging [1][2][3][4][5].Consequently, as shown in Fig. 5, our camera enables the possibility of performing imaging beyond the Abbe-Rayleigh criterion.In this case, we performed multiple experiments in which a superposition of partially distinguishable sources were imaged.The superposition was prepared using one coherent and one thermal light source.In Fig. 5a, we plot the predicted transverse separation s normalized by the Gaussian beam waist radius w 0 for both protocols.Here w 0 = λ/πNA, this parameter is directly obtained from our experiment.As demonstrated in Fig. 5a, our protocol enables one to resolve spatial features for sources with small separations even for diffraction-limited conditions.As expected for larger separation distances, the performance of our protocol matches the accuracy of intensity measurements.This is further demonstrated by the spatial profiles shown from Fig. 5b to d.The first row shows spatial profiles for three experimental points in Fig. 5a obtained through direct imaging whereas the images in the second row were obtained using our scheme for superresolving imaging.The spatial profiles in Fig. 5b show that both imaging techniques lead to comparable resolutions and the correct identification of the centroids of the two sources.However, as shown in Fig. 5c and d, our camera outperforms direct imaging when the separations decrease.Here, the actual separation is smaller than w 0 /2 for both cases.It is worth noticing that in this case, direct imaging cannot resolve spatial features of the sources.Here, the predictions of direct imaging become unstable and erratic.Remarkably, our simulations show an excellent agreement with the experimental data obtained for our scheme for superresolving imaging (see Methods section).
In conclusion, we demonstrated a robust quantum camera that enables superresolving imaging beyond the Abbe-Rayleigh resolution limit.Our scheme for quantum statistical imaging exploits the self-learning features of artificial intelligence to identify the statistical fluctuations of truly unknown mixtures of light sources.This particular feature of our scheme relies on a universal model based on the theory of quantum coherence to describe the photon statistics produced by the scattering of an arbitrary number of light sources.We demonstrated that the measurement of the quantum statistical fluctuations of photons enables one to overcome inherent limitations of existing superresolution protocols based on spatial mode projections [14,[18][19][20][21][22].We believe that our work represents a new paradigm in the field of optical imaging with important implications for microscopy, remote sensing, and astronomy [5][6][7][8][9][10][11].

Training of NN
For the sake of simplicity, we split the functionality of our neural network into two phases: the training and testing phase.In the first phase, the training data is fed to the network multiple times to optimize the synaptic weights through a scaled conjugate gradient backpropagation algorithm [42].This optimization seeks to minimize the Kullback-Leibler divergence distance between predicted and the real target classes [43,44].At this point, the training is stopped if the loss function does not decrease within 1000 epochs [45].In the test phase, we assess the performance of the algorithm by introducing an unknown set of data during the training process.For both phases, we prepare a data-set consisting of one thousand experimental measurements of photon statistics for each of the five classes.This process is formalized by considering different numbers of data points: 100, 500, ..., 9500, 10000.Following a standardized ratio for statistical learning, we divide our data into training (70%), validation (15%), and testing (15%) sets [46].The networks were trained using the neural network toolbox in MATLAB, which runs on a computer Intel Core i7-4710MQ CPU (@2.50GHz) with 32GB of RAM.

Fittings
To determine the optimal fits for Fig. 4e and f we design a search space based on Eqs. 13 and 3. To do so we first found the mean photon number of the input pixel, which will later be applied to constrain the search space.From here we allowed for the existence of up to three distinguishable modes which will be combined according to Eq. 3.Each of the modes contains an indistinguishable combination of up to one coherent and two thermal sources whose number distribution is given by Eq. 13.The total combination results in partially distinguishable combination and provides the theoretical model for our experiment.From here our search space is n=0 (p exp (n) − p th (n| n 1,t , n 2,t , n c )) 2 , where n i,t and n c are the mean photon numbers of that each thermal or coherent source contributes to each distinguishable mode respectively.The mean photon numbers of each source must add up to the experimental mean photon number, constraining the search.A linear search was then performed over the predicted mean photon numbers and the minimum was returned, providing the optimal fit.

Monte-Carlo Simulation of the Experiment
To demonstrate a consistent improvement over traditional methods, we also simulated the experiment using two beams, a thermal and a coherent, with Gaussian point spread functions over a 128×128 grid of pixels.At each pixel, the mean photon number for each source is provided by the Gaussian point spread function, which is then used to create the appropriate distinguishable probability distribution as given in Eq. 3, creating a 128×128 grid of photon number distributions.The associated class data for these distributions will then be fitted using to a set of pre-labeled disks using a genetic algorithm.This recreates our method in the limits of perfect classification.Each of these distributions is then used to simulate photon-number resolving detection.This data is then used to create a normalized intensity for the classical fit.We fit the image to a combination of Gaussian PSFs.This process is repeated ten times for each separation in order to average out fluctuations in the fitting.When combining the results of the intensity fits they are first divided into two sets.One set has the majority of fits return a single Gaussian, while the other returned two Gaussian the majority of the time.The set identified as only containing a single Gaussian is then set at the Abbe-Rayleigh diffraction limit, while the remaining data is used in a linear fit.This causes the sharp transition between the two sets of data.

FIG. 3 .
FIG.3.Projection of the feature space on the plane defined by the probabilities p(0), p(1), and p(2).The red points correspond to the photon statistics for coherent light, and the blue points indicate the photon statistics for thermal light fields.Furthermore, the brown dots represent the photon statistics produced by the scattering of two thermal light sources, and the black points show the photon statistics for a mixture of photons emitted by one coherent and one thermal source.The corresponding statistics for a mixture of one coherent and two thermal sources are indicated in green.As shown in a, the distributions associated to the multiple sources obtained for 10 data points are confined to a small region of the feature space.A similar situation prevails in b for 100 data points.As shown in panel c, the distributions produced with 1000 data points occupy different regions, although brown and black points keep closely intertwined.Finally, the separated distributions obtained with 10000 data points in d enable efficient identification of light sources.