Quantum limits of localisation microscopy

TL;DR

This paper establishes quantum limits on the precision of localizing multiple incoherent point sources in microscopy, revealing fundamental bounds and comparing measurement techniques against these limits.

Contribution

It introduces a quantum framework for localization microscopy, deriving bounds on estimation precision and comparing measurement strategies.

Findings

Quantum Fisher information bounds the localization precision.

Maximum independent parameters estimable is two, regardless of sources.

Spatial-mode demultiplexing can approach quantum limits.

Abstract

We show that localisation microscopy of multiple weak, incoherent point sources with possibly different intensities in one spatial dimension is equivalent to estimating the amplitudes of a classical mixture of coherent states of a simple harmonic oscillator. This enables us to bound the multi-parameter covariance matrix for an unbiased estimator for the locations in terms of the quantum Fisher information matrix, which we obtained analytically. In the regime of arbitrarily small separations we find it to be no more than rank two -- implying that no more than two independent parameters can be estimated irrespective of the number of point sources. We use the eigenvalues of the classical and quantum Fisher information matrices to compare the performance of spatial-mode demultiplexing and direct imaging in localisation microscopy with respect to the quantum limits.

Equations

$${\rho }_{\mathrm{opt}}\approx (1-ϵ){\rho }_{\mathrm{vac}}+ϵ\rho ,$$

$${\hat{c}}_i^{†}=\int dx{\Psi }_{\mathrm{PSF}}(x-{\chi }_i){\hat{c}}_x^{†},$$

$$\rho =\sum _{i=1}^N{w}_i|{\psi }_i⟩⟨{\psi }_i|,$$

$$|{\psi }_i⟩=\int dx{\Psi }_{\mathrm{PSF}}(x-{\chi }_i)|x⟩,$$

$${\psi }_{\mathrm{PSF}}(x)=\frac{1}{(2\pi {\sigma }^2{)}^{1/4}}{e}^{-\frac{x^2}{4{\sigma }^2}},$$

$$|{\psi }_i⟩=\sum _{k=0}^{\infty }\frac{{\alpha }_i^k}{\sqrt{k!}}{e}^{-{\alpha }_i^2/2}|{\varphi }_k⟩\equiv |{\alpha }_i⟩,$$

$$\rho \equiv {\rho }_{\mathbf{\alpha }}=\sum _{i=1}^N{w}_i|{\alpha }_i⟩⟨{\alpha }_i|,$$

$$\mathrm{Cov}[\mathbf{\alpha }]=\sum _{z}p(z|\mathbf{\alpha })(\mathbf{\alpha }-\tilde{\mathbf{\alpha }}{)}^T(\mathbf{\alpha }-\tilde{\mathbf{\alpha }}),$$

$$\mathrm{Cov}[\mathbf{\alpha }]\ge \frac{1}{Mϵ}{\left[\mathcal{C}({\rho }_{\mathbf{\alpha }},{\mathrm{\Pi }}_z)\right]}^{-1},$$

$$[\mathcal{C}({\rho }_{\mathbf{\alpha }},{\mathrm{\Pi }}_z){]}_{\mu \nu }=\sum _z\frac{1}{p(z|\mathbf{\alpha })}\frac{\partial p(z|\mathbf{\alpha })}{\partial {\alpha }_{\mu }}\frac{\partial p(z|\mathbf{\alpha })}{\partial {\alpha }_{\nu }}.$$

$$\underset{\{{\Pi }_z\}}{\max }\text{leavevmode}\mathcal{C}({\rho }_{\mathbf{\alpha }},{\mathrm{\Pi }}_z)\le \mathcal{Q}({\rho }_{\mathbf{\alpha }}),$$

$$[\mathcal{Q}({\rho }_{\mathbf{\alpha }}){]}_{\mu \nu }=\mathrm{Tr}\left[{\rho }_{\mathbf{\alpha }}\frac{L^{\mu }{L}^{\nu }+L^{\nu }{L}^{\mu }}{2}\right],$$

$$2\frac{\partial {\rho }_{\mathbf{\alpha }}}{\partial {\alpha }_{\mu }}=({\rho }_{\mathbf{\alpha }}{L}^{\mu }+L^{\mu }{\rho }_{\mathbf{\alpha }}).$$

$$\mathrm{Cov}[\mathbf{\alpha }]\ge \frac{1}{Mϵ}{\left[\mathcal{C}({\rho }_{\mathbf{\alpha }},{\mathrm{\Pi }}_z)\right]}^{-1}\ge \frac{1}{Mϵ}{\left[\mathcal{Q}({\rho }_{\mathbf{\alpha }})\right]}^{-1}.$$

$$\mathrm{Tr}\left({\rho }_{\mathbf{\alpha }}[L^{\mu },L^{\nu }]\right)=2\mathrm{Tr}\left(\mathrm{Im}\left({\rho }_{\mathbf{\alpha }}{L}^{\mu }{L}^{\nu }\right)\right)=0.$$

$${\mathcal{Q}}^{\prime }=B\mathcal{Q}{B}^T,$$

$${\rho }_{\mathbf{\alpha }}=A\left(\begin{array}{cc}{D}_{\mathbf{w}}& 0\\ 0& 0\end{array}\right)A^{†}\equiv A{\rho }_A{A}^{†}$$

$$A=\left(\begin{array}{ccccccc}|{\alpha }_1⟩& |{\alpha }_2⟩& \cdots & |{\alpha }_N⟩& {\hat{a}}^{†}|{\alpha }_1⟩& \cdots & {\hat{a}}^{†}|{\alpha }_N⟩\end{array}\right)$$

$${\rm Y}=A^{†}A,$$

$${\rm Y}=\left(\begin{array}{cc}{{\rm Y}}_{\alpha \alpha }& {{\rm Y}}_{\alpha d}\\ {{\rm Y}}_{d\alpha }& {{\rm Y}}_{dd}\end{array}\right),$$

$${\left({{\rm Y}}_{\alpha \alpha }\right)}_{ij}$$

$${\left({{\rm Y}}_{\alpha d}\right)}_{ij}={{\rm Y}}_{d\alpha }^{†}$$

$${\left({{\rm Y}}_{dd}\right)}_{ij}$$

$$=D_{\mathbf{\alpha }}{{\rm Y}}_{\alpha \alpha }{D}_{\mathbf{\alpha }}+{{\rm Y}}_{\alpha \alpha },$$

$${\partial }_j{\rho }_{\mathbf{\alpha }}=A{w}_j\left(\begin{array}{cc}-2{\alpha }_j{E}_j& E_j\\ E_j& 0\end{array}\right)A^T\equiv A({\partial }_j\rho {)}_A{A}^{†},$$

$$L^j=A{L}_A^j{A}^T=A\left(\begin{array}{cc}{L}_{\alpha \alpha }^j& L_{\alpha d}^j\\ L_{d\alpha }^j& L_{dd}^j\end{array}\right)A^T$$

$$2{\partial }_j{\rho }_A={\rho }_A{\rm Y}{L}_A^j+L_A^j{\rm Y}{\rho }_A,$$

$$[\mathcal{Q}({\rho }_{\mathbf{\alpha }}){]}_{jk}=\mathrm{Tr}\left({\partial }_j\rho L^k\right)=\mathrm{Tr}\left({\partial }_j{\rho }_A{\rm Y}{L}_A^k{\rm Y}\right).$$

$$\mathrm{vecb}\left(L_A^j\right)=\left[\begin{array}{c}|L_{\alpha \alpha }^j)\\ |L_{d\alpha }^j)\\ |L_{\alpha d}^j)\\ |L_{dd}^j)\end{array}\right],$$

$$2\mathrm{vecb}\left({\partial }_j{\rho }_A\right)=\left(\mathbb{I}\odot ({\rho }_A{\rm Y})+({\rho }_A{\rm Y})\odot \mathbb{I}\right)\mathrm{vecb}\left(L_A^j\right)$$

$$\mathrm{Tr}\left(A^{T}BC{D}^T\right)=\mathrm{vecb}{\left(A^T\right)}^T(D\odot B)\mathrm{vecb}\left(C\right),$$

$$=\mathrm{vecb}{\left({\partial }_i{\rho }_A\right)}^T({\rm Y}\odot {\rm Y})\mathrm{vecb}\left(L_A^j\right)$$

Full text

I Introduction

Precisely locating multiple single emitters is a key challenge in fluorescence microscopy. The process of estimating these locations depends on the quality of the image obtained by the microscope. One of the major limitations to the image quality, known since Abbe and Rayleigh, lies in spatially resolving objects substantially smaller than half the wavelength of the light involved Born and Wolf (1999). Known as the Rayleigh limit or diffraction limit, it is a consequence of the diffraction of light due to its wave nature.

Over the last couple of decades, ways to circumvent the Rayleigh limit in far-field fluorescence microscopy have been invented Thorley et al. (2014). Confocal methods such as STED, RESOLFT, and SSIM Hell and Wichmann (1994); Heintzmann et al. (2002); Gustafsson (2005); Hofmann et al. (2005) use patterned illumination to spatially modulate the fluorescence pattern of emitters within a diffraction-limited region such that not all of them emit simultaneously, thereby achieving sub-Rayleigh resolution. Other far-field methods such as PALM, fPALM and STORM Rust et al. (2006); Betzig et al. (2006); Hess et al. (2006) temporally modulate the fluorescence pattern of emitters with weak laser pulses stochastically such that only a low density of emitters are active within the Rayleigh limit at one time. Repeating the process many times, images with sub-Rayleigh resolution are reconstructed from the measured positions of individual emitters. These techniques, with resolution of tens of nanometers, have provided insights into biological processes at the cellular scale that were hitherto unattainable Huang et al. (2010).

Though immensely powerful and impressive, none of these methods seek to extract all the information available in the emitted light field. As in conventional fluorescence microscopy these techniques use ‘direct imaging’—intensity measurements on the image plane—to extract information from the incident light. That there is indeed more information in the light field to be extracted was shown by Tsang et al. (2016). Using methods from classical and quantum estimation theory, it was shown theoretically that two arbitrarily close incoherent point sources may be resolved, and that this may be achieved in practice using a spatial-mode demultiplexing (SPADE) measurement. In the few years since, theoretical studies have considered different source arrangements or parameters of interest Nair and Tsang (2016a); Lupo and Pirandola (2016); Kerviche et al. (2017); Chrostowski et al. (2017); Řehaček et al. (2017); Dutton et al. (2019) in one as well as in two and three spatial dimensions Ang et al. (2017); Yu and Prasad (2018); Napoli et al. (2019); Backlund et al. (2018). Other theoretical studies have explored various detection systems that could achieve the ultimate precision in imaging or get close to it Nair and Tsang (2016b); Yang et al. (2017); Rehacek et al. (2017); Řeháček et al. (2018). Several experiments have demonstrated some of the principles underlying these detection systems Tang et al. (2016); Paúr et al. (2016); Yang et al. (2016); Donohue et al. (2018); Paúr et al. (2018); Parniak et al. (2018); Zhou et al. (2019); Bonsma-Fisher et al. (2019). Advances in this area have been recently reviewed by Tsang (2019a).

Realistic imaging scenarios typically involve more than two point sources or even extended objects. It has been shown that an extended one-dimensional object much smaller than the Rayleigh limit described only in terms of its centroid and effective radius can be approximated by a two-level quantum system Chrostowski et al. (2017). Theoretical optimality of certain measurement techniques in estimating this effective radius size has also been established in one and two spatial dimensions Tsang (2017, 2018); Dutton et al. (2019). Order-of-magnitude bounds on the precision of estimating the normalised moments of extended sources smaller the Rayleigh limit have also been obtained Zhou and Jiang (2019); Tsang (2019b).

In this paper, we provide an analytical lower bound on an unbiased estimator’s covariance (mean square error) matrix for localisation microscopy – simultaneously estimating the locations of $N$ incoherent, weak point sources of unequal but known intensities in one spatial dimension. The bound is provided by the the quantum Fisher information matrix. For a Gaussian point spread function (PSF), we first describe the light field on the image plane as a classical mixture of coherent states. We use this to derive the quantum Fisher information matrix analytically. In the limit of the point sources approaching a single point, we find its rank to be no more than two. As the inverse of the quantum Fisher information matrix lower bounds the covariance matrix, our result implies that no more than two independent parameters can be estimated in localisation microscopy in the limit of arbitrarily small separations. In this limit, we provide a mathematical explanation for our observation in terms of an approximation of the light field involving only the first two Hermite-Gauss modes. Finally, we compare performance of conventional direct imaging and the recently proposed SPADE Tsang et al. (2016) in localisation microscopy with the quantum bounds we obtain. In the limit of the point sources approaching a single point, we find the classical Fisher information matrices for both these detection systems to be rank one. Furthermore, in the sub-Rayleigh limit, SPADE does not attain the quantum limit for localisation microscopy. For the subset of parameters where scalings may be optimal, we find SPADE to be short of the quantum limit in absolute precision.

This paper is organised as follows: In Section II we provide a quantum mechanical description of localisation microscopy. The appropriate framework to study the quantum limits of the localisation problem is quantum estimation theory, the toolbox of which is described in Section III, leading to the definition of the quantum Fisher information matrix (QFIM). In Section IV we provide an analytic expression of the QFIM for localisation microscopy, our main technical result. We then draw conclusions about its rank and its implications for localisation microscopy. We end in Section V with further insights and discussions about the sinc PSF and the potential of detection systems attaining the quantum limits of localisation microscopy.

II Quantum description of localisation microscopy

We consider localisation microscopy – the problem of estimating the locations of $N$ incoherent point sources or emitters located in a one-dimensional spatial configuration as in Fig. 1. As we assume them to be weak, such that on average no photons arrive on the image place within a coherence time with probability $(1-\epsilon)$, where $\epsilon\ll 1$ and one photon arrives with probability $\epsilon$. We also assume the optical field on the image plane to be quasi-monochromatic and paraxial Tsang et al. (2016). The quantum state of this optical field is then

[TABLE]

where we have neglected terms of second and higher orders in $\epsilon$ and $\rho_{\mathrm{vac}}=\ket{\mathrm{vac}}\bra{\mathrm{vac}}$ is the vacuum state and $\rho$ is the one-photon state.

The one-photon density matrix on the object plane is an incoherent mixture of position eigenstates $\rho=\sum_{i=1}^{N}w_{i}\ket{\chi_{i}}\bra{\chi_{i}}$, where $w_{i}$ are the relative intensities with $\sum_{i=1}^{N}w_{i}=1$. An imaging system maps $\hat{c}_{x}^{\dagger},$ the creation operator producing one photon in the position $x$ on the object plane, to the corresponding image plane operator $\hat{c}_{i}^{\dagger}$ Lupo and Pirandola (2016)

[TABLE]

where $\chi_{i}$ is the position on the source on the object plane and $\psi_{\mathrm{PSF}}(x)$ is the PSF. On the image plane this becomes

[TABLE]

where

[TABLE]

as follows from Eq. (2).

An ideal imaging system with $\Psi_{\mathrm{PSF}}(x)=\delta(x)$ is free of any Rayleigh limit as it transmits all spatial frequencies from the object to the image plane. In practice, a Gaussian PSF

[TABLE]

with $\sigma=\lambda/(2\pi\mathrm{NA}),$ where $\mathrm{NA}$ is the numerical aperture of the imaging system is a good approximation for quasimonochromatic paraxial light Zhang et al. (2006); Tsang et al. (2016) and also allows us to obtain analytical results. For such a PSF, the state of Eq. (3) has an intensity distribution of the form illustrated in Fig. 1. For a Gaussian PSF, the $\ket{\psi_{i}}$ can be expanded in the Hermite-Gauss (HG) basis as (See Appendix A)

[TABLE]

where $\ket{\phi_{k}}$ are the HG modes111Unlike the conventional quantum optical coherent states which reside in the phase space of the electromagnetic field, our coherent states reside in physical space on the image plane. This mathematical form was also identified by Dutton et al. (2019) but only used for numerical calculations.

This has the same mathematical form as the coherent states, produced by the displacement operator $\mathcal{D}(\alpha_{i})=e^{\alpha_{i}\hat{a}^{\dagger}+\alpha_{i}^{*}\hat{a}}$ Kok and Lovett (2010) acting on the ground state of the harmonic oscillator with $\alpha_{i}=\chi_{i}/2\sigma\in\mathbb{R}$ the dimensionless positions of the sources. Thus the one-photon state on the image plane is

[TABLE]

a classical mixture of coherent states in the HG basis.

The above is a quantum optical rendition of localisation microscopy—a classical optics problem. It enables us to harness the mathematical formalism associated with coherent states and provides a basis that spans the space of the quantum state as well as its derivative. The latter is an essential ingredient of deriving the quantum Fisher information matrix analytically in Section IV.1. We also hope that this description will provide insights into the quantum limits to localisation microscopy in the presence of shot noise and assist in designing detection systems that attain these quantum limits.

III Quantum Estimation Theory

Localisation has long been treated as an estimation problem with the unknown locations of the sources $\bm{\chi}\equiv\{\chi_{i}\},\leavevmode\nobreak\ i=1,\dots,N$ being the parameters to be estimated Ober et al. (2004); Chao et al. (2016). In our formulation, the limits to the localisation of the point sources are the same as estimating the amplitudes $\bm{\alpha}\equiv\{\alpha_{i}\},\leavevmode\nobreak\ i=1,\dots,N$ of the coherent states in Eq. (7). Let these estimates be $\tilde{\bm{\alpha}}\equiv\{\tilde{\alpha}_{i}\}$. The precision of our estimate is then given by the covariance (or mean square error) matrix defined as

[TABLE]

where $p(z|\bm{\alpha})$ is the probability distribution of the collected data labelled by, for instance, the pixel $z$ on the image plane. $\operatorname{Cov}[\bm{\alpha}]$ is a positive symmetric matrix whose $i-$th diagonal element denotes the variance of an estimator of $\alpha_{i}$ given the data collected. The $(i,j)$-th off-diagonal element denote the covariance in the estimation of $\alpha_{i}$ and $\alpha_{j}$.

Given the data collected, the maximum amount of information that can be extracted from it to obtain the most precise estimate of the locations is given by the Cramér-Rao bound (CRB) Cover and Thomas (2006). For an unbiased estimator, this bound is given by

[TABLE]

where $M$ is the number of coherence times over which the data is collected, making $M\epsilon$ the total photon count. This inequality is saturable but generally only in the asymptotical limit of many repetitions Trees (2001). $\mathcal{C}(\rho_{\bm{\alpha}},\mathrm{\Pi}_{z})$ is the classical Fisher information matrix (CFIM) whose elements are given by Paris (2009)

[TABLE]

The probability distribution $p(z|\bm{\alpha})=\mathrm{Tr}\left(\rho_{\bm{\alpha}}\mathrm{\Pi}_{z}\right)$ results from detecting the light on the image plane using a specific detection system $\mathrm{\Pi}_{z}$. Fluorescence microscopy typically employs intensity detectors \mathrm{\Pi}_{z}=\{\mbox{|n\rangle!\langle n|}_{z}\},n=0,1,\cdots, at each pixel $x,$ known as direct imaging. It is then evident that the CFIM depends on the detection system used, and not surprising that it determines the amount of information that can be extracted from the light field at the image plane.

To identify the quantum limit on the precision of localisation microscopy, the CFIM must be maximised over all possible physically allowed detection systems. This set is given by positive operator-valued measures (POVMs) Nielsen and Chuang (2000) and the maximisation is bounded as Holevo (1982); Helstrom (2009).

[TABLE]

by the quantum Fisher information matrix (QFIM). Its matrix elements are given by

[TABLE]

with $L_{\mu}$ being the symmetric logarithmic derivative (SLD) corresponding to the parameter $\alpha_{\mu}$. The SLD is determined by the Lyapunov equation

[TABLE]

The quantum limit to localisation microscopy is thus given by

[TABLE]

The QFIM depends only on the light field on the image plane and determines the maximum amount of information that can be extracted from it using detection systems allowed by quantum mechanics. Deriving an analytical expression for $\mathcal{Q}(\rho_{\bm{\alpha}})$ for state in Eq. (7) is our main result, which we present in the next section.

A practical issue following the identification of the quantum limit is its attainability. For cases where a single parameter is unknown then a measurement can be found to satisfy the equality of Eq. (11), which involves projecting onto the eigenstates of the SLD Braunstein and Caves (1994); Paris (2009). However this strategy does not generalise to multiple parameters, as in localisation microscopy, in general.

For multi-parameter estimation the attainability is tantamount to saturating the second inequality in Eq. (14). A necessary condition for the saturability of any scalar form of Eq. (14) is the satisfaction of weak commutativity Matsumoto (2002); Ragy et al. (2016)

[TABLE]

Moreover, through the quantum theory of asymptotic normality Kahn and Guţă (2009), this condition becomes sufficient with the application of collective measurements over multiple copies of $\rho_{\bm{\alpha}}$  Matsumoto (2002); Ragy et al. (2016).

Any scalar function of the covariances can be bounded by the inverse of QFIM with the lower bound following from the spectral decomposition of QFIM. To that end, calculating the eigenvalues of the QFIM and their scaling is of importance for the multi-parameter estimation. For localisation microscopy, $\rho_{\bm{\alpha}}$ as in Eq. (7) as well as its derivative are real matrices. Thereby, $L^{\mu},L^{\nu}$ are also real and the above condition is always satisfied222Since the density matrix and its derivatives are real-valued in the orthonormal $\{\ket{\phi_{k}}\}$ basis, Eq. (13) is a system of equations with real coefficients. Hence $L^{\nu}$ must be real as well, and so $\rho_{\bm{\alpha}}L^{\mu}L^{\nu}$ is real-valued. We thank Ben Wang for bringing this to our attention. The quantum limit for localisation microscopy is therefore attainable, at least in principle, although collective measurements over multiple copies Matsumoto (2002); Ragy et al. (2016) of the light field on the image plane may be required.

Alternative parameterisations of the system—where the new parameters $\bm{\alpha}^{\prime}$ are functions of the old parameters $\bm{\alpha}$—can be dealt with by a transformation of the QFIM. Given the transformation matrix $B$ with elements $B_{ij}=\partial\alpha_{i}/\partial\alpha^{\prime}_{j}$, the QFIM of the transformed parameters is Paris (2009)

[TABLE]

provided the transformation is non-singular. This can be used to recast our results in terms of, for instance, the moments of the point source distribution.

IV Results

We now present our main result – the analytical expression of the QFIM for localisation microscopy. This expression allows us to conclude that the QFIM is a rank two matrix as $\alpha_{i}\rightarrow 0$. Eq. (14) then implies that the eigenvalues of $\operatorname{Cov}[\bm{\alpha}]$ remains finite for no more than two independent parameters. Thus, no more than two independent parameters can be estimated from the entire set $\bm{\alpha}$ as $\bm{\alpha}\rightarrow 0$.

We lack a fully satisfactory physical explanation for this restriction on the number of estimable parameters, but provide an explanation involving only the first two Hermite-Gauss modes for $\alpha_{i}\ll 1$.

IV.1 Analytical expression of QFIM

The state in Eq. (7) can be expressed in the basis of $\{\ket{\alpha_{i}},\hat{a}^{\dagger}\ket{\alpha_{i}}\}$ as

[TABLE]

where

[TABLE]

and $D_{\bm{w}}=\operatorname{diag}\left(w_{1},w_{2},\cdots,w_{N}\right)$ denotes a diagonal matrix. Although the basis used in Eq. (17) is non-orthogonal this representation can still be used to evaluate the QFIM Genoni and Tufarelli (2019). The coherent states $\{\ket{\alpha_{i}}\}$ are linearly independent and span the support of the state in Eq. (7). The support of the derivative is spanned by $\{\ket{\alpha_{i}}\}$ and $\{\hat{a}^{\dagger}\ket{\alpha_{i}}\}$, which are also linearly independent.

The Grammian matrix

[TABLE]

whose elements consist of the scalar products between the basis vectors $\braket{\alpha_{j}}{\alpha_{k}}$, $\braket{\alpha_{j}}{\hat{a}^{\dagger}}{\alpha_{k}}$, $\braket{\alpha_{j}}{\hat{a}}{\alpha_{k}}$, and $\braket{\alpha_{j}}{\hat{a}\hat{a}^{\dagger}}{\alpha_{k}}$ is in block form,

[TABLE]

where

[TABLE]

and $D_{\bm{\alpha}}=\operatorname{diag}\left(\alpha_{1},\alpha_{2},\cdots,\alpha_{n}\right).$

Since $\partial_{\alpha}\ket{\alpha}=(\hat{a}^{\dagger}-\alpha)\ket{\alpha}$ for real $\alpha$, the derivative of the quantum state is

[TABLE]

where $\partial_{j}$ denotes the derivative with respect to $\alpha_{j}$ and $(E_{j})_{kl}=\delta_{jk}\delta_{jl}$. Similarly, the SLD $L_{A}^{i}$ can be written in the generic form

[TABLE]

where $L_{\alpha\alpha}^{j}$ corresponds to the elements $\bra{\alpha_{i}}L^{j}\ket{\alpha_{j}}$, $L_{\alpha d}^{j}$ to $\bra{\alpha_{j}}L^{j}\hat{a}^{\dagger}\ket{\alpha_{i}}$ etc. The Lyapunov equation Eq. (13) can be now rewritten as

[TABLE]

and the QFIM elements from Eq. (12) as

[TABLE]

Using the Tracy-Singh block kronecker product $\odot$ and the block column "$\mathrm{vecb}$" operator Koning et al. (1991) defined as

[TABLE]

where $|X)=\mathrm{vec}(X)$ is the column vectorisation of a matrix and $(X|$ its transpose. Eq. (24) can be blockwise vectorised to

[TABLE]

with $\mathbb{I}$ being the identity matrix. Using the matrix identity Koning et al. (1991)

[TABLE]

the QFIM elements from Eq. (25) can be re-expressed as

[TABLE]

where we have defined

[TABLE]

which is the outstanding quantity to be determined.

We now recast Eq.(27) and (29) as

[TABLE]

Putting it all together, we obtain

[TABLE]

where

[TABLE]

and $\{\upsilon_{\alpha\alpha},\upsilon_{\alpha d},\upsilon_{d\alpha},\upsilon_{dd}\}$ defines the inverse of $\Upsilon$ via

[TABLE]

Note that the inverse $\Upsilon^{-1}$ always exists since $\Upsilon$ is the Grammian matrix of linearly independent vectors. The elements of $\Upsilon^{-1}$ can be found using the formula of blockwise inversion (See Appendix B).

Noticing that $|\Gamma_{dd}^{j})$ does not contribute in Eq. (29), Eq. (32) can be reduced to

[TABLE]

where $\mathbb{A}$ is a $3N^{2}\times 3N^{2}$ invertible matrix unless $\alpha_{i}=\alpha_{j}$ for some $i,j$, which is a singular case for which the rank of the density matrix reduces. Hence the unique solution to Eq. (35) is

[TABLE]

where the block matrices that compose the $\mathbb{A}^{-1}$ can be found by using the formulas for blockwise inversion (See Appendix B).

Substituting Eq. (36) into Eq. (29) gives us the QFIM elements

[TABLE]

where $S^{-1}=\left(\Upsilon_{\alpha\alpha}^{-1}\otimes D_{\bm{w}}+D_{\bm{w}}\otimes\Upsilon_{\alpha\alpha}^{-1}\right)^{-1}$ is an $N^{2}\times N^{2}$ matrix and $\Upsilon_{\alpha\alpha}$ is the inverse of the submatrix of $\Upsilon$ which exists, as it is the Grammian matrix of linearly independent vectors $\{\ket{\alpha_{i}}\}$. Eq. (37) is an analytic expression for the QFIM elements for localisation microscopy and our main result.

Fig. 2 shows the elements of the QFIM for the localisation microscopy of three point sources. We choose them to be equidistant, that is, $(\alpha_{1},\alpha_{2},\alpha_{3})=(x,2x,3x)$ and $w_{1}=w_{2}=w_{3}=1/3$ for illustration purposes. Note the non-zero off-diagonal elements evidencing correlations in the precision around and below the Rayleigh limit of $x\sim 1.$

While the diagonal elements are all non-vanishing, more crucially as $x\rightarrow 0$ the diagonal and off-diagonal elements combine to make the QFIM singular. This is revealed by a closer analysis of the QFIM matrix as in Fig. 3 which shows that only two of its eigenvalues remain non-zero as the sources approach each other. This is in spite of all the diagonals elements of the QFIM remaining non-zero even as $x\rightarrow 0$, as Fig. 2 shows.

This behaviour of only two non-zero eigenvalues also holds for other values of $N$. We have explicitly checked this for $N=4,\dots,10$ as well as when the sources are not equally spaced. In Fig. 7 in Appendix B we plot the eigenvalues of the QFIM for $N=4,5$ as further examples. In the case of different relative intensities the results are the same except of the limiting case of one extremely bright source $w_{j}\gg 1,w_{i\neq j}\ll 1$, where the rank of the QFIM is approximately one (Fig. 8 in Appendix B).

Since the QFIM has rank two as $x\rightarrow 0,$ its inverse is ill-defined except on a two-dimensional subspace. This implies that the $N\times N$ covariance matrix for localisation microscopy, as per Eq. (14), will also be unbounded except on a two-dimensional subspace. Thus, no more than two independent parameters can be estimated in localisation microscopy as the point sources approach each other.

In other words, the rank-deficient nature of the QFIM shows that a form of the Rayleigh limit resurfaces for any $N>2$. This had been suggested by previous works based on order-of-magnitude bounds for the diagonal elements on the CFIM Zhou and Jiang (2019) or uppers bounds on the diagonal elements of the QFIM Tsang (2019b). Our analytical expression for the full QFIM—its diagonal and off-diagonal elements for any $N$—shows that this rank two behaviour is truly quantum mechanical in origin. Furthermore, knowing the full QFIM matrix allows us to uncover the nature in which $N-2$ of the eigenvalues approach zero. We return to the behaviour in which this rank deficiency or Rayleigh limit emerges in Sec. V.

IV.2 Why rank two?

We now provide an explanation for the rank deficiency of the QFIM in the regime of small separations which can be seen as the re-emergence of the Rayleigh limit. To that end, we expresses the state in Eq. (7) in terms of the real-valued displacement operator $\mathcal{D}(\alpha_{i})=e^{\alpha_{i}\hat{c}_{i}^{\dagger}-\alpha_{i}\hat{c}_{i}}$ as

[TABLE]

In the limit of very small separations ($\alpha_{i}\ll 1$), the displacements are approximately

[TABLE]

where $\mathbb{I}$ is the identity operator and the displacement $\alpha_{i}$ is real. Up to the second order in $\alpha_{i}$, the normalised quantum state of the light field on the image place is then

[TABLE]

where $\mathcal{C}_{i}$ are the first two moments

[TABLE]

Eq. (40) describes the state of two-level quantum system—the two levels being the first two HG modes. A similar approximation which described the state relative to a PSF centred at a fixed reference point was used in Ref. Chrostowski et al. (2017) to estimate the centroid and the effective radius of a distribution of incoherent point sources. We now consider the more general problem of estimating the location of $N$ point sources.

The QFIM for $\bm{\alpha}$ (See Appendix C) is

[TABLE]

with

[TABLE]

where $M_{11}=(\mathcal{C}_{2}-1)\mathcal{C}_{2}$, $M_{12}=M_{21}=\mathcal{C}_{1}(1-2\mathcal{C}_{2})$, $M_{22}=4\mathcal{C}_{1}^{2}-1$, $\mathcal{A}=\left(\mathcal{C}_{2}-1\right)\mathcal{C}_{2}+\mathcal{C}_{1}^{2}$, and $\bm{I}=\left(1\;1\dots 1\right)^{T}$.

The QFIM $\mathscr{Q}$ is an $N\times N$ matrix, which is a product of three matrices of dimensions $N\times 2$, $2\times 2$ and $2\times N$. Since $\operatorname{rank}(AB)\leq\min\left\{\operatorname{rank}(A),\operatorname{rank}(B)\right\}$, and the matrix $\mathcal{M}$ has rank 2, the QFIM $\mathscr{Q}$ has rank no more than two. Although a two-level quantum system has the potential of estimating three real parameters, localisation microscopy in this limit can estimate only two as the two-level system possesses a real density matrix333As the localisation parameters $\bm{\alpha}$ are real, $\mathrm{Tr}\left(\rho_{\bm{\alpha}}^{(2)}\,\sigma_{y}\right)=2\,\mathrm{Im}\left(\mathcal{C}_{1}\right)=0,$ where $\sigma_{y}$ is the Pauli $Y$ matrix. This is another way of arguing that as the point sources get closer, the light field on the image plane has enough information to estimate only two parameters. A physical reason for this observation would be highly desirable.

V Discussion

Our analytical expression for the QFIM for localisation microscopy has enabled us to show that as point sources get closer, no more than two independent parameters can be estimated. A rank-deficient QFIM occurs when the quantum state does not contain enough information to permit the estimation some of the parameters or combinations thereof. The parameters that can be estimated correspond to the non-zero eigenvalues of the QFIM. Without additional knowledge of the source distribution this restricts us to estimating functions of the first two moments $f(\mathcal{C}_{1},\mathcal{C}_{2})$ only deep in the sub-Rayleigh limit. As Eq. (42) shows, when all $\{\alpha_{j}\}$ are unknown as in localisation microscopy, there is vanishing information about any single $\alpha_{i}$ itself. This is in contrast to the scalar QFI $[\mathcal{Q}(\rho_{\bm{\alpha}})]_{ii}$ for $\alpha_{i}$ which is non zero, but assumes that all the other $\{\alpha_{j}\}$ are known. The manner in which the eigenvalues of the QFIM tend to zero is of interest in the search for optimal detection systems for localisation microscopy. Numerical fitting in Fig. 5 shows the vanishing eigenvalues of the QFIM approach zero polynomially. The degree of the polynomial is given by $d=2\lfloor\frac{\mu-1}{2}\rfloor$, where $\mu$ is the order the eigenvalue when arranged in descending order and $\lfloor\cdot\rfloor$ is the floor function. These scalings are now extracted from the elements of the full QFIM of the localisation parameters $\bm{\alpha}$ – rather than from bounds on estimating the various moments independently as in previous works Zhou and Jiang (2019); Tsang (2019b).

Unlike the latter, we can now compare the absolute performance of detection systems for localisation microscopy relative to its quantum limit. Indeed, while Fig. 5 shows the $2n$-th eigenvalue of the QFIM closely parallel to the $n$-th eigenvalue of the CFIM for SPADE Tsang et al. (2016), there is a large gap in the absolute terms. This could be due to the sub-optimality of SPADE for estimating the $\lfloor N/2\rfloor$ parameters it is sensitive to444Conventional SPADE is not sensitive to all the parameters needed to describe the sources’ distribution, only its even moments Tsang (2019b); Zhou and Jiang (2019); Tsang (2019a).. Similar scalings were observed with detection using superpositions of the conventional SPADE basis Tsang (2017, 2018); Zhou and Jiang (2019) that are sensitive to the other half of the moments. For reference over a range of separations, Fig. 10 in Appendix E shows the eigenvalues of the CFIM for SPADE as well as direct imaging. Note that for both, the CFIM tends towards a rank one matrix.

Finally, although our analytical result is derived with a Gaussian PSF, we expect the rank deficiency of the QFIM to be present in a more general family of PSFs. To that end, Fig. 6 shows the numerically obtained eigenvalues of the QFIM for three equidistant point sources of equal intensities under a sinc PSF (See Appendix D) defined as

[TABLE]

This PSF is the exact form for diffraction through a sharp one-dimensional slit which in its principal peak is well-approximated as Gaussian.

An approximation involving the first two spherical Bessel modes as in Sec. IV.2 can be performed for a sinc PSF as well, leading to similar insights. A proof of this rank deficiency for arbitrary PSFs and a physical explanation remains an open question.

To conclude, we have obtained several insights into the quantum limits of localisation microscopy via an analytical expression for the QFIM. In particular, the behaviour of the eigenvalues of the QFIM deep in the sub-Rayleigh limit revealed that only two parameters are eventually estimable. It also enabled us to compare the performance of known detection systems relative to the quantum limit in absolute terms, a question left open in the literature Tsang (2019a). The gap identified by us should motivate the search for detection systems, ideally on a single copy of the light field on the image plane, seeking to reduce or eliminate it.

VI Acknowledgements

We thank Francesco Albarelli, Jamie Friel, Mankei Tsang, Liang Jiang, Lijian Zhang, Alex Retzker for enlightening discussions. This study has been supported by the UK EPSRC (EP/K04057X/2), the UK National Quantum Technologies Programme (EP/M01326X/1, EP/M013243/1) and the University of Warwick Global Partnership Fund.

Appendix A Expressing the density matrix in the HG basis

The density matrix is written in terms of the kets $\ket{\psi_{i}}$, which are expressed in the position space as in Eq. (3). We assume a normalised Gaussian point spread function (PSF) of the form

[TABLE]

and so

[TABLE]

The kets $\ket{\psi_{i}}$ can be expressed in terms of the complete Hermite-Gauss modes as

[TABLE]

where $\ket{\phi_{q}}$ are the Hermite-Gauss modes, which can be expressed in the position space as Tsang et al. (2016)

[TABLE]

where $\mathrm{H}_{q}(x)$ are the Hermite polynomials. The coefficients of the expansion Eq. (A.3) are

[TABLE]

Setting $\frac{\chi_{i}}{2\sigma}=\alpha_{i}$ we get

[TABLE]

which has the same mathematical form as the coherent states with $\{\ket{\phi_{q}}\}$ forming the Fock basis Kok and Lovett (2010).

The state in Eq. (3) can be also written in terms of the displacement operators $\mathcal{D}(\alpha_{i})=e^{\alpha_{i}(a^{\dagger}-a)}$, with $\alpha_{i}=\frac{\chi_{i}}{2\sigma}\in\mathbb{R}$

[TABLE]

where $\mathcal{D}(\alpha)$ is the displacement operator.

The derivative of each coherent state with respect to its real amplitude $\alpha$ is given by

[TABLE]

which yields the formula Eq. (22).

Appendix B Analytic results for $N$ sources

The Tracy-Singh product Tracy and Singh (1972); Koning et al. (1991) defined for matrices $A$ and $B$ subdivided into blocks $A_{ij}$ and $B_{kl}$ is $A\odot B$ where the $(i,j)$-th block of $A\odot B$ is $A_{ij}\odot B$ whose $(k,l)$-th block is in turn $A_{ij}\otimes B_{kl}$. That is if $A,B$ are block matrices with

[TABLE]

then the Tracy-Singh product is

[TABLE]

Using the above definition, the matrix of Eq. (31) is found to be

[TABLE]

where the elements of $\Upsilon^{-1}$ can be found using the formula of blockwise inversion:

[TABLE]

with

[TABLE]

The inverse of the block matrix $\Upsilon_{\alpha\alpha}$ exists, because it is the Gramian matrix of the linear independent vectors $\ket{\alpha_{i}}$.

For the QFIM elements we need to evaluate the inverse of the top left $3N^{2}\times 3N^{2}$ part of the matrix of Eq. (B.1) which we denote $\mathbb{A}$. In order to obtain the inverse of $\mathbb{A}$, we need to further partition $\mathbb{A}$ as

[TABLE]

with

[TABLE]

The inverse of $\varpi$ is

[TABLE]

The elements of $\mathbb{A}^{-1}$ will be given by the formulas

[TABLE]

After calculations and by substituting the $\Upsilon^{-1}$ elements from Eq. (B.3), we derive the explicit form of $\mathbb{A}^{-1}$ elements:

[TABLE]

The QFIM elements are then obtained from Eq. (36) and (B.8)

[TABLE]

Finally, to complement the discussion in the main text,we present some further examples of the QFIM eigenvalues for $N=4,5$ spurces and in Fig. 8 we present the eigenvalues of the QFIM for 3 sources in the case of unequal weights (relative intensities) Fig.(8).

Appendix C Analytic results for $x_{i}\ll\sigma$

The state in the sub-diffraction regime is given by Eq. (40). The derivative can be calculated immediately from this formula and it is

[TABLE]

By solving the SLD equation $\partial_{a_{i}}\rho_{\bm{\alpha}}=(\rho_{\bm{\alpha}}L_{i}+L_{i}\rho_{\bm{\alpha}})$, we can determine the SLDs in the $\{\ket{0},\ket{1}\}$ basis:

[TABLE]

Knowing the SLDs, we can obtain the QFIM of Eq. (42).

As already mentioned in the main text, the rank of the QFIM only depends on the matrix

[TABLE]

of Eq. (42), with the elements of this matrix given by Eq. (43). The eigenvalues $\mu_{1},\mu_{2}$ of the matrix Eq. (C.3) are

[TABLE]

The condition for the eigenvalues to be zero is

[TABLE]

The first part $0=3\mathcal{C}_{1}^{2}-1-\sqrt{\left(3\mathcal{C}_{1}^{2}-1\right)^{2}}$ is always true, as it reduces to the identity $\left(3\mathcal{C}_{1}^{2}-1\right)^{2}=\left(3\mathcal{C}_{1}^{2}-1\right)^{2}$. For the second part we have

[TABLE]

Substituting $\mathcal{C}_{2}$ and $\mathcal{C}_{1}$ we get

[TABLE]

Since $\alpha_{i}$ are strictly positive, except one that can be zero, this sum of positive terms cannot be equal to zero. Therefore, this statement is always false. Thus, the Eq. (C.5) becomes $(1\land 0)\lor(1\land 0)=0$, which means that the two eigenvalues can never be zero and the QFIM will be rank 2.

Appendix D Calculation of the QFI for the Sinc PSF

The expansion of the Sinc function on the HG modes is not ideal for numerical calculations. Instead we use the spherical Bessel function of the $1^{\mathrm{st}}$ kind and express the states onto those modes in which we then truncate. If the PSF is a $\operatorname{sinc}$ function, the $\ket{\psi_{i}}$ are

[TABLE]

We can use the identity Abramowitz and Stegun (1964)

[TABLE]

where $J_{q}(x)$ is the spherical Bessel function of the $1^{st}$ kind. The spherical Bessel function are orthogonal in all $\mathbb{R}$

[TABLE]

therefore we can define the orthonormal basis

[TABLE]

The set of the spherical Bessel functions is a basis in $\mathbb{R}$, but is not complete since it is not a resolution of identity as we can see from Eq. (D.2). Hence, we can expand the sinc function on the bessel function basis, using the identity Eq. (D.2):

[TABLE]

Using the identity for the Bessel functions

[TABLE]

we can also have an expression for the derivative of $\ket{\psi_{i}}$

[TABLE]

We see that both the state $\rho$ and its derivatives are completely expressed within the basis $\ket{j_{q}}$. This means that we can use the definition of the SLD (Eq. D.8) and express the SLD in the same basis.

[TABLE]

In this way the fact that the specific basis is not complete does not affect our calculations.

For the numerical calculations we have to truncate our state in the appropriate amount of modes. From Figs. 6 and 9, we can see that our conclusions do not change with the use of a non-Gaussian PSF.

Appendix E Eigenvalues of the CFIM for SPADE and Direct Imaging

Finally, we present the eigenvalues of the CFIM for SPADE and direct imaging fir a large range of separations.