Accuracy of the Gaussian Point Spread Function model in 2D localization microscopy

Sjoerd Stallinga; Bernd Rieger

doi:10.1364/OE.18.024461

1. Introduction

The exploitation of the physics of time dependent fluorescence has enabled the emergence of the field of far-field optical nanoscopy [1, 2]. This time dependence can be achieved by on/off switching of fluorophores, which is either intrinsically stochastic (blinking) or deliberately engineereed (photoswitching, saturation or depletion). Within the spectrum of optical nanoscopy techniques [3–9] localization microscopy holds the promise of ease-of-use because it uses a conventional widefield fluorescence microscopy set-up. One of the key factors in the success of localization microscopy is the Point Spread Function (PSF) model used for fitting the emitter positions in the sparse images, as the achievable localization accuracy is affected by the correctness of the PSF model. So far, mainly Gaussian model PSFs have been used, because of their conceptual simplicity and computational efficiency, despite the lack of physical foundation. High objective numerical aperture values require full vectorial modelling of the PSF. In addition, the orientation of the radiating dipoles plays a significant role in the observable spot shape. In many cases the emitter can either rotate or change its conformation freely during the lifetime of the excited state of the emitter, so for many excitation-emission cycles an average over randomly distributed emission dipole orientations will be observed. In other cases we must consider ensembles of fixed but randomly oriented dipoles, like fluorescent beads, and we need to average over the different dipole orientations as well, but now weighted with the excitation efficiency as a function of dipole orientation. In yet a third class of cases, the single emitter cannot rotate freely, such as in solid-state or cryofreezed samples or for emitters immobilized at a surface or in a matrix, and the fixed dipole orientation needs to be taken into account in the analysis, as it affects the observed spot shape [10, 11].

The question how to accurately model the PSF in accordance with the laws of optics has been addressed by a number of authors. The starting point is the theory laid down by Wilson, Higdon, Török, and co-workers [12–15] and adapted to single molecule localization by Enderlein and co-workers [16] and recently by Mortensen and co-workers [17]. The latter two references restrict their analysis to the case of 2D Total Internal Reflection Fluorescence (TIRF) imaging without any aberrations.

The goal of this paper is to find out how good or bad the popular Gaussian PSFs really are. We will look at the effects of dipole orientation and of the conventional optical aberrations (defocus, astigmatism, coma, spherical aberration). To that end the model of [16,17] is generalized to include arbitrary aberrations. This model is used for generating ‘ground truth’ simulated images, which are subsequently used as input for a Maximum Likelihood Estimator (MLE) based on a Gaussian model PSF. In ref. [18] and in a recent paper co-authored by one of us [19] it has been shown that such an estimator achieves the Cramér-Rao Lower Bound (CRLB) in the presence of background and Poissonian noise, down to photon counts well below 100. Here we analyze how the Gaussian-based MLE deviates from the theoretical CRLB due to model errors. For the sake of simplicity we will also restrict our attention to 2D imaging. This entails the study of only the lateral localization accuracy, the different methods for axial localization are not taken into account in the model and the analysis. It does not entail the absence of defocus, i.e. we do consider deviations of the imaged 2D-plane from the focal plane. Our analysis is therefore appropriate to TIRF-imaging but also applies to conventional widefield imaging.

2. Vectorial theory of imaging of a dipole emitter

2.1. Preliminaries and definitions

The emitter dipole is located in a medium with refractive index n_med adjacent to a cover slip with refractive index n_cov. The emitter is imaged with an immersion objective lens with numerical aperture NA_ob designed for an immersion fluid with refractive index n_imm. The objective lens is assumed to be corrected for focusing onto the interface between the cover slip and the medium. The intersection of the optical axis with this interface is taken to be the origin of the coordinate system in object space. The emitted radiation is collected by the objective lens with focal length F_ob and focused by the tube lens with focal length F_im onto the detector. Both lenses are assumed to be aplanatic and the imaging system is assumed to be telecentric, i.e. the aperture stop is located at the back focal plane of the objective lens, which coincides with the front focal plane of the tube lens. The magnification of the imaging system is M = F_im/F_ob, thus making the numerical aperture in image space NA_im = NA_ob/M. In practice M ≫ 1, so that NA_im ≪ 1. The radius of the stop is given by R = F_obNA_ob = F_imNA_im. Throughout this paper we use scaled coordinates. The pupil coordinates are scaled with the stop radius R scaling the pupil to the unit circle, the object and image coordiates are scaled with the diffraction lengths λ/NA_ob (object space) and λ/NA_im (image space). We use ‘u’ for object and image coordinates and ‘v’ for Fourier/pupil coordinates throughout this text. The emitter dipole vector is oriented along the unit-vector d⃗ = (sin θ_d cos ϕ_d, sin θ_d sin ϕ_d, cos θ_d), where θ_d is the polar angle, and ϕ_d is the azimuthal angle, and has magnitude d₀. The (2D) position of the emitter with respect to the focal point is r⃗_d = (x_d, y_d). The scaled position is u⃗_d = NA_obr⃗_d/λ.

2.2. Model for the PSF for arbitrary aberrations

The electric field in the pupil plane at point v⃗ = (v_x, v_y, 0) corresponds to the plane wave in object space with wavevector along (sin θ_med cos ϕ, sin θ_med sin ϕ, cos θ_med), so:

v_{x} = \frac{n_{med} sin θ_{med}}{N A_{ob}} cos ϕ,

v_{y} = \frac{n_{med} sin θ_{med}}{N A_{ob}} sin ϕ .

The electric field in the pupil plane has only non-zero values for the Cartesian components j = x,y given by [12–15]:

E_{pub, j} (\vec{v}, {\vec{u}}_{d}) = \frac{π d_{0}}{ɛ_{0} λ^{2} F_{ob} {(1 - v^{2} {NA}_{ob}^{2} / n_{med}^{2})}^{1 / 4}} \sum_{k = x, y, z} q_{j k} (\vec{v}) d_{k} e^{2 π i {\vec{u}}_{d} \cdot \vec{v}},

and with the polarization vectors:

q_{x k} = cos ϕ T_{p} p_{k} - sin ϕ T_{s} s_{k},

q_{y k} = sin ϕ T_{p} p_{k} + cos ϕ T_{s} s_{k},

which depend on the p and s basis polarization vectors p⃗ = (cos θ_med cos ϕ, cos θ_med sin ϕ, –sin θ_med) and s⃗ = (–sin ϕ, cos ϕ, 0) and on the Fresnel coefficients:

T_{a} = T_{a, med - cov} T_{a, cov - imm},

for a = p, s and where the Fresnel coefficients for the two contributing interfaces (see in particular [14] for details of the effects of stratified media) are defined by:

T_{a, 1 - 2} = \frac{2 c_{a, 1}}{c_{a + 1} + c_{a + 2}},

for a = p,s and with c_p,l = n_l/ cosθ_l and c_s,l = n_l cos θ_l for l = med, cov, imm.

The field in the pupil plane is modified by the pupil function

C (\vec{v}) = {\begin{matrix} A (\vec{v}) exp (i W (\vec{v})), & | \vec{v} | \leq 1, \\ 0, & | \vec{v} | > 1, \end{matrix}}

with A(v⃗) the amplitude and W(v⃗) the phase (aberration) of the pupil function, and the field in the image plane (on the detector) at scaled position u⃗ is found by subsequently applying a Fourier transform. As the numerical aperture of the tube lens is typically rather small compared to unity we may safely apply scalar diffraction theory to the x and y components of the field seperately. This gives:

E_{im, j} (\vec{u}, {\vec{u}}_{d}) = E_{0} \sum_{k = x, y, z} w_{j k} (\vec{u} - {\vec{u}}_{d}) d_{k},

with:

E_{0} = \frac{π^{2} {NA}_{ob} {NA}_{im} d_{0}}{i ɛ_{0} λ^{3}},

and:

w_{j k} (\vec{u}) = \frac{1}{π} \int d^{2} v C (\vec{v}) \frac{q_{j k} (\vec{v})}{{(1 - v^{2} {NA}_{ob}^{2} / n_{med}^{2})}^{1 / 4}} exp (- 2 π i \vec{u} \cdot \vec{v}) .

The measured signal (PSF) is proportional to the z-component of the Poynting-vector. In view of the low numerical aperture on the image side this is simply proportional to the square of the modulus of the electric field. The PSF then follows as:

I (\vec{u}, {\vec{u}}_{d}) = \frac{1}{2} ɛ_{0} c \sum_{j = x, y} {| E_{im, j} (\vec{u}, {\vec{u}}_{d}) |}^{2} = I_{0} \sum_{k, l = x, y, z} [\sum_{j = x, y} w_{j k} (\vec{u} - {\vec{u}}_{d}) w_{j l} {(\vec{u} - {\vec{u}}_{d})}^{*}] d_{k} d_{l} .

with:

I_{0} = \frac{π^{4} c {NA}_{ob}^{2} {NA}_{im}^{2} {| d_{0} |}^{2}}{2 ɛ_{0} λ^{6}} .

If there are no aberrations present the equations of refs. [16, 17] are reproduced as a special case. The current theory generalizes on these prior results by including arbitrary aberrations, in particular asymmetric ones like astigmatism and coma.

The total power captured by a pixel centered at u⃗_p is found by integrating over the pixel area A_p:

P ({\vec{u}}_{p}, {\vec{u}}_{d}) = \frac{λ^{2}}{{NA}_{im}^{2}} \int_{A_{p}} d^{2} u I (\vec{u}, {\vec{u}}_{d}) = P_{0} \sum_{k, l = x, y, z} [\sum_{j = x, y} \int_{A_{p}} d^{2} u w_{j k} (\vec{u} - {\vec{u}}_{d}) w_{j l} {(\vec{u} - {\vec{u}}_{d})}^{*}] d_{k} d_{l},

with:

P_{0} = \frac{π^{4} c {NA}_{ob}^{2} {| d_{0} |}^{2}}{2 ɛ_{0} λ^{4}} .

The dipole magnitude d₀ is given by:

d_{0} = α {\vec{E}}_{ill} \cdot {\vec{d}}_{ex},

with α the polarizability, E⃗_ill the electric field of the illumination and d⃗_ex a unit-vector along the excitation dipole. For a fixed dipole the excitation and emission dipole orientation are equal: d⃗_ex = d⃗. So, apart from the bilinear dependence on dipole orientation explicit in Eqs. (12) and (14) the PSF also depends on dipole orientation implicitly via the excitation efficiency.

There is an interesting relation between the theory presented here and the theory of Richards and Wolf and many subsequent authors for the vectorial field in the focal region when a polarized plane wave in the pupil plane is focused by the objective lens. The exact same functions w_jk(u⃗) turn out to determine the field close to focus. The electric field in the focal plane can be written as (see [20, 21] and references therein):

E_{foc, k} (\vec{u}) = \frac{π F_{ob}}{i λ} \sum_{j = x, y} A_{j} w_{j k} (\vec{u}),

where (A_x, A_y) is the unit-vector representing the polarization of the plane wave in the pupil plane.

2.3. Orientational averaging

There are two cases in which we must average over dipole orientation. The first case deals with fluorophores that have a rotational diffusion time that is much smaller than the lifetime of the excited state. In that case the emission dipole is fully uncorrelated to the excitation dipole and an independent average over both the excitation and emission dipole orientations must be taken. The second case in which averaging over dipole orientation is needed is the case of an ensemble of fixed but randomly oriented fluorescent molecules, for example as in a fluorescent bead. Now the excitation and emission dipoles are the same and the excitation probability poportional to |E⃗_ill · d⃗|² must be taken into account as a weight factor in the averaging prccedure over the different dipole orientations contributing to the PSF. In general this will lead to a different expression for the PSF. This has been explicitly shown for the TIRF-case by Mortensen et al. [17]. In the appendix we derive explicit expressions for the different orientationally averaged PSFs, for arbitrary states of polarization and coherence of the illumination electric field, and for arbitrary aberrations of the optical system. In the following we only report on results for fixed dipoles and for temporal orientationally averaged dipoles, the case of ensemble orientationally averaged dipoles is not considered any further from hereon.

The fixed and averaged dipole cases described here can be seen as limiting cases of a general distribution over dipole orientations characterized by a Probability Distribution Function (PDF), as described in ref. [22]. The fixed dipole case corresponds to an infinitely sharply peaked PDF, the averaged dipole case to a uniform PDF. Generalizing to an arbitary PDF would allow for a description of the effects of e.g. an arbitrary rotational diffusion time or of thermally driven fluctuations around a preferential orientation.

2.4. The case with azimuthal symmetry

In case of azimuthally symmetric aberrations (defocus, spherical aberration) the functions w_jk(u⃗) can be expressed as:

w_{x x} (u, ψ) = F_{0} (u) + F_{2} (u) cos (2 ψ),

w_{x y} (u, ψ) = F_{2} (u) cos (2 ψ),

w_{x z} (u, ψ) = i F_{1} (u) cos ψ,

w_{y x} (u, ψ) = F_{2} (u) cos (2 ψ),

w_{y y} (u, ψ) = F_{0} (u) - F_{2} (u) cos (2 ψ),

w_{y z} (u, ψ) = i F_{1} (u) sin ψ,

with the three integrals over functions containing the Bessel functions J_k(x), k = 0, 1, 2 well-known from the theory of Richards and Wolf and generalizations thereof [20]:

F_{0} (u) = 2 \int_{0}^{1} d v \frac{v (T_{s} + T_{p} \sqrt{1 - v^{2} {NA}_{ob}^{2} / n_{med}^{2}})}{2 {(1 - v^{2} {NA}_{ob}^{2} / n_{med}^{2})}^{1 / 4}} J_{0} (2 π u v) exp (i W (v)),

F_{1} (u) = 2 \int_{0}^{1} d v \frac{v^{2} T_{p} {NA}_{ob} / n_{med}}{{(1 - v^{2} {NA}_{ob}^{2} / n_{med}^{2})}^{1 / 4}} J_{1} (2 π u v) exp (i W (v)),

F_{2} (u) = 2 \int_{0}^{1} d v \frac{v (T_{s} - T_{p} \sqrt{1 - v^{2} {NA}_{ob}^{2} / n_{med}^{2}})}{2 {(1 - v^{2} {NA}_{ob}^{2} / n_{med}^{2})}^{1 / 4}} J_{2} (2 π u v) exp (i W (v)) .

The resulting electric field in the image plane turns out to be:

E_{im, x} (u, ψ) = E_{0} [F_{0} (u) sin θ_{d} cos ϕ_{d} + F_{2} (u) sin θ_{d} cos (2 ψ - ϕ_{d}) + i F_{1} (u) cos θ_{d} cos ψ,]

E_{im, y} (u, ψ) = E_{0} [F_{0} (u) sin θ_{d} sin ϕ_{d} + F_{2} (u) sin θ_{d} sin (2 ψ - ϕ_{d}) + i F_{1} (u) cos θ_{d} sin ψ],

and the resulting intensity is:

\begin{array}{l} I (u, ψ) & = & I_{0} [({| F_{0} (u) |}^{2} + {| F_{2} (u) |}^{2}) {sin}^{2} θ_{d} + {| F_{1} (u) |}^{2} {cos}^{2} θ_{d} \\ + 2 Im {(F_{0} (u) + F_{2} (u)) F_{1} {(u)}^{*}} sin (2 θ_{d}) cos (ψ - ϕ_{d}) \\ - 2 Re {F_{0} (u) F_{2} {(u)}^{*}} {sin}^{2} θ_{d} cos (2 ψ - 2 ϕ_{d})], \end{array}

In the absence of aberrations, this expression for the PSF is equivalent to the results of [16, 17]. For a dynamic dipole orientation the different terms average out to give:

I (u) = \frac{1}{3} I_{0} (2 {| F_{0} (u) |}^{2} + 2 {| F_{2} (u) |}^{2} + {| F_{1} (u) |}^{2}),

which does not depend on the azimuthal angle ψ anymore.

3. The effect of aberrations and dipole orientation on spot shape and localization accuracy

It is worthwile to study the expressions for the PSF for the fixed and free dipole cases in some detail before discussing the results of the numerical simulations.

The PSF for a freely rotating dipole is centrosymmetric but deviates significantly from the classical Airy-shape of scalar diffraction theory for the high NA values used in practice. Actually, the PSF effectively approaches a Gaussian shape with small shoulders at the location of the first Airy-rings (see Fig. 1). If the PSF is considered in the window of width 2 λ / NA it may be approximated with even better accuracy by a Gaussian plus a background [17,19]. This change in PSF-shape is due to the various vectorial diffraction effects that smear out the fringe structure of the Airy-distribution. For a water immersion objective the best-fit Gaussian has a standard deviation approximately equal to 0.25 λ /NA, which is equal to one pixel for Nyquist sampling. The approximation of the PSF with a Gaussian has been analyzed in greater detail in ref. [23].

Fig. 1 Dipole orientation averaged PSF for a 100 ×/1.25 water immersion objective, an Airy-distribution for the same NA, and a Gaussian with standard deviation σ = 0.25 λ /NA.

Download Full Size | PDF

For the case at hand, determining the location of a single emitter it may be anticipated that the use of a Gaussian model PSF does not introduce a very harmful model error. Even for symmetric aberrations like defocus and spherical aberration no large loss of localization accuracy is expected. However, for asymmetric aberrations like coma we may reasonably expect a significant impact on localization accuracy.

The PSF for a fixed emission dipole is generally non-centrosymmetric [16, 17]. The spot shape can change dramatically of character for different polar dipole angles. For a planar dipole (θ_d close to π/2) the first and fourth term in Eq. (29) dominate. The first term is centrosymmetric, the fourth term makes the spot slighlty elliptical, similar to image spots suffering from astigmatism. For an axial dipole (θ_d close to zero) the second term dominates, giving rise to a doughnut spot (dark spot in the middle with bright ring), similar to a spot with significant defocus and spherical aberration. For intermediate polar angles the third term gains importance. This term makes the spot asymmetric and this resembles the effect of coma, when the spot becomes asymmetric with an increased side lobe in one direction.

Figure 2 shows simulated image spots for different fixed dipole orientations and for different aberrations. The qualitative similarity discussed before between the spot distortions that occur for different polar dipole angles and the spot distortions that occur for different types of aberrations is quite apparent. It is also clear that the similarity is only qualitative, as the spot shapes cannot be matched very well. The relation between dipole orientation induced spot distortions and conventional aberrations such as defocus, astigmatism, coma, and spherical aberration suggests that aberrations may be used to cancel the dipole orientation effects. Alas, it turns out that this does not work very well; aberrations tend to make things worse in all circumstances.

Fig. 2 Simulated spot shapes for different dipole orientations and aberrations. Left: Image spots for NA_ob = 1.45 (which is larger than n_med = 1.33) for a fixed dipole with polar angle θ_d = 0 (a), θ_d = π/4 (b), θ_d = π /2 (c). Right: Image spots for a free dipole aberrated by 300 m λ RMS defocus (d), 160 m λ RMS coma (e), and 100 m λ RMS astigmatism (f).

Download Full Size | PDF

It may be anticipated that similar conclusions may be drawn regarding the symmetry or asymmetry of the spot, and type of distortions of spot shape, for more general optical imperfections, such as the polarization aberrations described by McGuire and Chipman [24, 25].

The main source of error in emitter localization can be expected to come from the asymmetry term. It appears from Eq. (29) that this term only contributes if the function F₁(u) and/or the function F₀(u) + F₂(u) is complex. This can occur when the Fresnel coefficients are complex. This happens for a numerical aperture larger than the medium refractive index, e.g. when an oil immersion objective is used to collect light from emitters in a watery medium. Then evanescent waves can couple from the dipole emitter into the objective, and for these evanescent waves the Fresnel-coefficients are complex. This is exactly the case studied in [16, 17], and indeed they have noticed the strong asymmetry due to the third term in the PSF. This dipole orientation induced asymmetry in the PSF can be eliminated by stopping down the objective to a numerical aperture smaller than the medium refractive index at the expense of a reduction in photon collection efficiency. The photon collection efficiency increases with NA_ob because the solid angle subtended by the aperture of the microscope objective increases, and for NA_ob > n_med also by evanescent wave coupling into the marginal region of the objective aperture. This so-called supercritical fluorescence hugely increases the photon output [16] and the ratio of photon collection efficiency between NA_ob = 1.45 and NA_ob = 1.25 turns out to be a factor of about 3. It follows that either the acquisition time is increased with this factor of 3 or the statistical localization error is increased with a factor $\sqrt{3}$ . This seems a small price to pay for avoiding large systematic localization errors due to spot shape distortions. Naturally, the amplitude of the evanescent wave contributions decay exponentially with the distance from the the cover slip. As a consequence, the spot asymmetry effect and the supercritical fluorescence based photon output gain only apply to the TIRF context, and vanish for emitters more than a few wavelengths away from the cover slip.

The PSF can become markedly asymmetric, even for the case with NA_ob < n_med. This effect, observed by Bartko and Dickson [10], occurs for non-zero aberrations, in particular defocus. Then the functions F_j(u) also become complex valued. Fig. 3 shows the combined effect of a fixed dipole orientation and a defocus of 72 m λ RMS (which corresponds to Marechal’s diffraction limit). For λ = 500 nm, NA_ob = 1.25, and n_med = 1.33 it can be estimated with the methods of ref. [26] that this amount of defocus corresponds to a change in focal position in the medium of about 150 nm. The asymmetry effect is quite pronounced with a shift in the spot peak of about 0.16 λ /NA_ob ≈ 64 nm. The typical error is found from this worst case error by multiplying with the average of the factor sin(2θ_d) over all appearing dipole orientations. This gives a factor 2/3 making the typical error 43 nm. This error is a systematic error that must be added to the statistical error related to the finite number of collected photons. Longer measurements in order to decrease the statistical error therefore cease to be meaningful beyond a certain photon count. From these considerations it may be anticipated that a fixed dipole orientation implies that the spot shape is more prone to distortions caused by aberrations. Clearly, this effect reduces aberration tolerances for 2D localization microscopy and complicates attempts to make 3D localization microscopy sufficiently accurate.

Fig. 3 Simulated spot shape for NA_ob = 1.25 (which is less than n_med = 1. 33) for a fixed dipole with polar angle θ_d = π /4, without aberrations (a) and with 72 m λ RMS defocus (b).

Download Full Size | PDF

4. Numerical results

We performed numerical simulations for analyzing the accuracy of using a Gaussian model PSF by applying the MLE described in [19] to ‘ground truth’ images calculated with the vector imaging theory outlined in the previous sections. The field and intensity in the image plane were calculated by evaluating the 2D Fourier transform in Eq. (11) using the chirp z-transform method [27] with a 118 × 118 grid for the pupil plane and a 11 × 11 grid for the image plane. The image plane grid corresponds to a block of 11 × 11 pixels, the pixel size was taken to be equal to λ /4NA_im (Nyquist sampling). For a water immersion objective with NA_ob = 1.25 and a wavelength λ = 500 nm this implies a pixel size in object space of 100 nm. It was checked that the simulated error did not depend significantly on pixel size in the range [0.15–0.30] × λ /NA_im. In addition, it appeared that convolution over the finite pixel size did not have a significant impact on the error either, so subsequently it was ignored in the simulation. Random configurations were generated with emitter positions distributed with a normal distribution centered on the central pixel of the block and with standard deviation of one pixel. The dipole orientation for each configuration was randomly selected from a uniform distribution over the full 4π solid angle. After calculation of the ‘ground truth’ PSF an image was generated by distributing N photons over the entire image following Poisson statistics with the ‘ground truth’ PSF as rate parameter. The MLE with Gaussian model PSF was subsequently used to find the emitter locations as well as the best fit to the background and to the height and width of the Gaussian. About 20 iterations were sufficient for convergence of all fitted parameters. In addition to simulations with fixed, randomly selected, dipole orientation we have also done simulations with freely rotating dipoles giving the rotationally averaged PSF. The statistical error was evaluated from 500 random configurations. For a small photon count (less than 20 in the 11 × 11 block of pixels) outliers with a fitted photon count more than twice the ‘real’ photon count were removed before calculating the localization error. We wish to emphasize that similar results to ours are to be expected from other estimators than the MLE given a sufficient photon count, i.e. a good signal-to-noise-ratio [28]. In particular, Least Mean Squares (LMS) methods for estimating the center of gravity of the spot [29], possibly in conjunction with spatial high-pass filtering for background removal, should be mentioned in this context.

4.1. Effect of dipole orientation

Figure 4 shows the effect of dipole orientation on the accuracy of the MLE-fit. The estimated error is plotted as a function of the detected number of photons for both fixed and free dipole orientation and for NA-values below and above the medium refractive index n_med = 1.33 (NA = 1.25 and 1.45). The CRLB that can be achieved with the MLE-estimator if the actual PSF was a Gaussian with a σ of one pixel is also plotted in the figure. It appears that a free dipole orientation gives rise to a localization accuracy that approaches the Gaussian based CRLB for all photon counts. This may be attributed to the centrosymmetric nature of the PSF and the fact that the MLE-estimator not only estimates the emitter location but also the PSF width and background. A fixed but unknown dipole orientation hardly changes the performance of the estimator provided the NA is smaller than the medium refractive index, as then the asymmetry term in the PSF is not present. If the NA is larger than the medium refractive index, however, the error of the Gaussian model PSF is so large that this will lead to a systematic localization error approaching 10 nm that must be added to the statistical error related to the finite photon count. It may be concluded that a Gaussian model PSF gives an excellent fit of emitter location provided the emitter dipoles can rotate freely during the excited state lifetime of the emitter or provided the objective NA is less than the medium refractive index in case the dipole orientation is fixed. In practice, the latter requirement boils down to using water immersion objectives instead of oil immersion objectives.

Fig. 4 Estimated error as a function of detected number of photons for fixed/free dipole orientation and NA_ob values below and above the refractive index of the medium, and the theoretical Gaussian based CRLB.

Download Full Size | PDF

4.2. Effect of aberrations

Figure 5 shows the error as a function of detected number of photons for free and fixed dipoles and for different aberration settings. We have investigated the effect of defocus, astigmatism, coma, and spherical aberration, taking values of 50% and 100% of Marechal’s diffraction limit (72 mλ RMS). An optical system is considered well-corrected if the total RMS aberation is less than the diffraction limit. In the simulations NA_ob = 1.25 and we take Nyquist sampling, so for the wavelength λ = 500 nm the pixel size is equal to the Nyquist unit λ /4NA_ob = 100 nm. For freely rotating dipoles aberrations hardly have an effect on localization accuracy, with the exception of coma, which gives rise to a systematic localization error of about several nm for the diffraction limit. For dipoles with a fixed orientation the situation is entirely different. Now the same level of aberrations introduce systematic localization errors (the plateau for large photon counts) of several tens of nm. For example, for a diffraction limited amount of defocus the x and y systematic error for diffraction limited defocus (the plateau for large photon count) is about 23 nm, so the total mean square error is 33 nm, in reasonable agreement with the estimate based on the average peak offset of the simulated spot shapes equal to 43 nm. Astigmatism and spherical aberration have a similar effect albeit with a smaller systematic localization error for the same level of aberrations. Naturally, coma gives foremost an effect in the asymmetry direction with systematic errors of 22 nm for the diffraction limited level.

Fig. 5 Estimated error as a function of detected number of photons for different aberration settings. Left (a to d): averaged over dipole orientation, right (e to f): randomly selected fixed dipole orientation. First row: defocus (a and e), second row: astigmatism (b and f), third row: coma (c and g), fourth row: spherical aberration (d and h).

Download Full Size | PDF

4.3. Effect of background

We have noticed several subtleties regarding the fitted background. For a small photon count the fitted background is zero, for high photon count the fitted background is essentially non-zero. Somewhat suprisingly, the transition between the two regimes does not appear to be continuous. For intermediate photon counts, around 100, the MLE fit routine can converge to two qualitatively different local minima, one with zero background and one with non-zero background. The relative occurence of the two types of solutions gradually changes from the zero background type to the non-zero background type as the photon count increases. It appears that there is no correlation between the type of solution that is found and the dipole orientation or any other parameter. The small shoulder in the intermediate regime of the curves of Fig. 5 that can develop when aberrations are present is attributed to this effect. For the case of coma, the asymmetric aberration, another complication arises. Now the zero background solution appears to have a systematic offset whereas the non-zero background solution has not. For that reason the offset averaged over all configurations varies depending on the relative occurence of the two types of solutions, i.e. depending on photon count. This will introduced another systematic error if the variations in photon count per emitter are large. These are largely determined by shot noise variations, although the effect of dipole orientation on detection efficiency may be quite significant for the fixed dipole case. For a TIRF setup it is therefore advisable to use a circularly polarized illumination, as this minimizes the dependence of excitation efficiency on dipole orientation.

Simulations with ‘ground truth’ images with an added background reveal another effect of interest. For an added background of one photon per pixel the bistable effect described above has disappeared, but now a new complication is observed. For small photon counts (total number of ‘signal’ photons in the block of 11 × 11 pixels is less than the total number of background photons) the MLE fit routine converges to Gaussians with a significantly larger width than the physical PSF width and a significant fraction of the background photons is counted as signal. We believe this is due to randomly distributed shot noise peaks in the background that are mistakenly identified as the image spot. As a consequence, the localization uncertainty is larger than the CRLB for a Gaussian PSF with the width of the physical PSF, even if the effect of the added background is incorporated in the CRLB-curve [19]. The increase in localization error in this regime is evident from Fig. 6. The deviations of the MLE curve from the CRLB in the low photon count regime can possibly be described by efficiency metrics such as the Barankin bound [30] or even more elaborate measures [31].

Fig. 6 Estimated error as a function of detected number of photons for fixed dipole orientation and NA_ob = 1.25, for different defocus settings in presence of a background of one photon per pixel, and the theoretical Gaussian based CRLB with this background.

Download Full Size | PDF

5. Conclusion

It may be concluded from the foregoing analysis that the use of a Gaussian PSF does not compromise the localization accuracy provided the emitter dipole orientations are not fixed, even for aberrations that are within the usual tolerances of optical design, with the exception of coma. For fixed emitter dipole orientations the use of a Gaussian PSF is no longer justified. In case the aberrations are corrected well below the diffraction limit the use of a water immersion objective instead of an oil immersion objective is quite beneficial, but for aberration levels that can occur in practice the image spots are distorted to such an extent that the use of Gaussian model PSFs introduces systematic errors of up to about 40 nm, even for water immersion objectives.

There are several ways forward to alleviate the problem for fixed dipoles. First of all, we may follow the lines of [17] and use a fully vectorial PSF model, including aberrations, in the estimator for the emitter position. Second, we can also try a computationally efficient set of basic functions to expand the PSF, and use the expansion coefficients to fit a possibly asymmetric spot shape. Possibly, a limited set of extra shape parameters can be accommodated, at the expense of a moderate increase in the photon count required for a good quality fit of all parameters involved. Finally, the experimental setup can be modified with polarization optics [10, 11, 15] for measuring dipole orientation parameters instead of estimating them.

Aberration induced spot deformations are also anticipated to be problematic for 3D localization microscopy, outside of the context of TIRF imaging. The measurement of the axial coordinate requires the introduction of astigmatism with a cylinder lens [32, 33] or making two images with different defocus [34, 35] (see also [36] for an experimental comparison). In addition, the different axial positions introduce large amounts of defocus. This level of aberrations introduces spot distortions that are so large that they must be dealt with in the position estimator with a more sophisticated model PSF than the simple, straightforward, and convenient Gaussian.

A. Temporal and ensemble orientational averaging

Here we show how to do the averaging for both temporal and ensemble orientational averaging, for arbitrary states of polarization and coherence of the illumination electric field, and for arbitrary aberrations of the optical imaging system. For the sake of notational convenience we drop all functional dependencies from the equations and write Eq. (12) as:

I = \sum_{i j k l} C_{i j} A_{k l} d_{ex, i} d_{ex, j} d_{k} d_{l},

with the matrix A_kl defined by:

A_{k l} = I_{0} \sum_{l} w_{j l} w_{k l}^{*},

and with the mutual coherence:

C_{i j} = {〈 Re (E_{ill, i} E_{ill, j}^{*}) 〉}_{E},

where the angular brackets 〈...〉_E denotes the averaging over the stochastic occurences of the electric field (state of coherence). The electric field of the illumination can always be written as:

{\vec{E}}_{ill} = \frac{1}{\sqrt{2}} A_{0} (\vec{p} + i \vec{q}),

with A₀ a (complex) amplitude, and p⃗ and q⃗ two real unit vectors that span the polarization ellipse. For example, for a linear polarization p⃗ and q⃗ are equal, while for a circular polarization they are mutually orthogonal. The mutual coherence matrix can be expressed in terms of these vectors as:

C_{i j} = \frac{1}{2} {| A_{0} |}^{2} {〈 p_{i} p_{j} + q_{i} q_{j} 〉}_{E} .

For fully incoherent illumination we have $C_{i j} = ({| A_{0} |}^{2} / 3) δ_{i j}$ .

For the temporal orientational averaging we need independent averaging over the excitation and emission dipole, which gives a PSF:

〈 I 〉 = \frac{1}{9} \sum_{k l} C_{k k} A_{l l},

where it is used that:

〈 d_{j} d_{k} 〉 = \frac{1}{3} δ_{j k},

with δ_jk Kronecker’s delta function (and a similar expression for the excitation dipole). For the ensemble orientational averaging we have d⃗_ex = d⃗ and we can use the orientational average over the product of four dipole factors d_k:

〈 d_{i} d_{j} d_{k} d_{l} 〉 = \frac{1}{15} (δ_{i j} δ_{k l} + δ_{i l} δ_{j k} + δ_{i k} δ_{j l}),

(this expression can be easily derived from symmetry arguments), in order to obtain the averaged intensity:

〈 I 〉 = \frac{1}{15} \sum_{k l} (C_{k k} A_{l l} + 2 C_{k l} A_{k l}) .

The second term between brackets makes this ensemble orientationally averaged PSF generally different in character from the temporal orientationally averaged PSF. For the special case of fully incoherent illumination, however, it is found that both cases of orientational averaging give the same PSF.

References and links

1. D. Evanko, “Primer: fluorescence imaging under the diffraction limit,” Nat. Methods 6, 19–20 (2009). [CrossRef]

2. S. W. Hell, “Microscopy and its focal switch,” Nat. Methods 6, 24–32 (2009). [CrossRef] [PubMed]

3. S. W. Hell, “Far-Field Optical Nanoscopy,” Science 316, 1153–1158 (2007). [CrossRef] [PubMed]

4. E. Betzig, G. H. Patterson, R. Sougrat, O. W. Lindwasser, S. Olenych, J. S. Bonifacino, M. W. Davidson, J. Lippincott-Schwartz, and H. F. Hess, “Imaging Intracellular Fluorescent Proteins at Nanometer Resolution,” Science 313, 1643–1645 (2006). [CrossRef]

5. K. A. Lidke, B. Rieger, T. M. Jovin, and R. Heintzmann, “Superresolution by localization of quantum dots using blinking statistics,” Opt. Express 13, 7052–7062 (2005). [CrossRef] [PubMed]

6. M. J. Rust, M. Bates, and X. Zhuang, “Sub-diffraction-limit imaging by stochastic optical reconstruction microscopy (STORM),” Nat. Methods 3, 793–795 (2006). [CrossRef] [PubMed]

7. J. Flling, M. Bossi, H. Bock, R. Medda, C. A. Wurm, B. Hein, S. Jakobs, C. Eggeling, and S. W. Hell, “Fluorescence nanoscopy by ground-state depletion and single-molecule return,” Nat. Methods 5, 943–945 (2008). [CrossRef]

8. A. Egner, C. Geisler, C. von Middendorff, H. Bock, D. Wenzel, R. Medda, M. Andresen, A. C. Stiel, S. Jakobs, C. Eggeling, A. Schnle, and S. W. Hell, “Fluorescence nanoscopy in whole cells by asynchronous localization of photoswitching emitters,” Biophys. J. 93, 3285–3290 (2007). [CrossRef] [PubMed]

9. M. Heilemann, S. van de Linde, M. Schttpelz, R. Kasper, B. Seefeldt, A. Mukherjee, P. Tinnefeld, and M. Sauer, “Subdiffraction-resolution fluorescence imaging with conventional fluorescent probes,” Angew. Chem. Int. Ed. Engl. 476172–6176 (2008). [CrossRef] [PubMed]

10. A. P. Bartko and R. M. Dickson, “Imaging three-dimensional single molecule orientations,” J. Phys. Chem. B 103, 11237–11241 (1999). [CrossRef]

11. P. Dedecker, B. Muls, J. Hofkens, J. Enderlein, and J. Hotta, “Orientational effects in the excitation and de-excitation of single molecules interacting with donut-mode laser beams,” Opt. Express 15, 3372–3383 (2007). [CrossRef] [PubMed]

12. T. Wilson, R. Juskaitis, and P. D. Higdon, “The imaging of dielectric point scatterers in conventional and confocal polarisation microscopes,” Opt. Commun. 141, 298–313 (1997). [CrossRef]

13. P. Torök, P. D. Higdon, and T. Wilson, “Theory for confocal and conventional microscopes imaging small dielectric scatterers,” J. Mod. Opt. 45, 1681–1698 (1998). [CrossRef]

14. O. Haeberlé, M. Ammar, H. Furukawa, K. Tenjimbayashi, and P. Török, “The point spread function of optical microscopes imaging through stratified media,” Opt. Express 11, 2964–2969 (2003). [CrossRef] [PubMed]

15. M. R. Foreman, C. M. Romero, and P. Török, “Determination of the three-dimensional orientation of single molecules,” Opt. Lett. 33, 1020–1022 (2008). [CrossRef] [PubMed]

16. J. Enderlein, E. Toprak, and P. R. Selvin, “Polarization effect on position accuracy of fluorophore localization,” Opt. Express 14, 8111 (2006). [CrossRef] [PubMed]

17. K. I. Mortensen, L. S. Churchman, J. A. Spudich, and H. Flyvbjerg, “Optimized localization analysis for single-molecule tracking and super-resolution microscopy,” Nat. Methods 7, 377–381 (2010). [CrossRef] [PubMed]

18. R. J. Ober, S. Ram, and E. S. Ward, “Localization Accuracy in Single-Molecule Microscopy,” Biophys. J. 86, 1185–1200 (2004). [CrossRef] [PubMed]

19. C. S. Smith, N. Joseph, B. Rieger, and K. A. Lidke, “Fast, single-molecule localization that achieves theoretically minimum uncertainty,” Nat. Methods 7, 373–375 (2010). [CrossRef] [PubMed]

20. S. Stallinga, “Axial birefringence and the light distribution close to focus,” J. Opt. Soc. Am. A 18, 2846–2859 (2001). [CrossRef]

21. S. Stallinga, “Light distribution close to focus in biaxially birefringent media,” J. Opt. Soc. Am. A 21, 1785–1798 (2004). [CrossRef]

22. M. R. Foreman, S. S. Sherif, and P. Török, “Photon statistics in single molecule orientational imaging,” Opt. Express 15, 13597–13606 (2007). [CrossRef] [PubMed]

23. B. Zhang, J. Zerubia, and J.-C. Olivo-Marin, “Gaussian approximations of fluorescence microscope point-spread function models,” Appl. Opt. 46, 1819–1829 (2007). [CrossRef] [PubMed]

24. J. P. McGuire Jr. and R. A. Chipman, “Polarization aberrations. 1. Rotationally symmetric optical systems,” Appl. Opt. 33, 5080–5100 (1994). [CrossRef] [PubMed]

25. J. P. McGuire Jr. and R. A. Chipman, “Polarization aberrations. 2. Tilted and decentered optical systems,” Appl. Opt. 33, 5101–5107 (1994). [CrossRef] [PubMed]

26. S. Stallinga, “Compact description of substrate-related aberrations in high numerical aperture optical disk readout,” Appl. Opt. 44, 949–958 (2005). [CrossRef]

27. J. L. Bakx, “Efficient computation of optical disk readout by use of the chirp z transform,” Appl. Opt. 41, 4879–4903 (2002). [CrossRef]

28. A. van den Bos, Parameter Estimation for Scientists and Engineers (Wiley & Sons, New Jersey, 2007). [CrossRef]

29. R. E. Thompson, D. R. Larson, and W. W. Webb, “Precise nanometer localization analysis for individual fluorescent probes,” Biophys. J. 82, 2775–2783 (2002). [CrossRef] [PubMed]

30. M. F. Kijewski, S. P. Müller, and S. C. Moore, “The Barankin bound: a model of detection with location uncertainty,” Proc. SPIE 1768, 153–160 (1992). [CrossRef]

31. S. P. Müller, C. K. Abbey, F. J. Rybicki, S. C. Moore, and M. F. Kijewski, “Measures of performance in nonlinear estimation tasks: prediction of estimation performance at low signal-to-noise ratio,” Phys. Med. Biol. 50, 3697–3715 (2005). [CrossRef] [PubMed]

32. L. Holtzer, T. Meckel, and T. Schmidt, “Nanometric three-dimensional tracking of individual quantum dots in cells,” Appl. Phys. Lett. 90, 053902 (2007). [CrossRef]

33. B. Huang, W. Wang, M. Bates, and X. Zhuang, “Three-dimensional super-resolution imaging by stochastic optical reconstruction microscopy,” Science 319, 810–813 (2008). [CrossRef] [PubMed]

34. E. Toprak, H. Balci, B. H. Blehm, and P. R. Selvin, “Three-dimensional particle tracking via bifocal imaging,” Nano Lett. 7, 2043–2045 (2007). [CrossRef] [PubMed]

35. M. F. Juette, T. J. Gould, M. D. Lessard, M. J. Mlodzianoski, B. S. Nagpure, B. T. Bennett, S. T. Hess, and J. Bewersdorf, “Three-dimensional sub-100nm resolution fluorescence microscopy of thick samples,” Nat. Methods 5, 527–530 (2008). [CrossRef] [PubMed]

36. M. J. Mlodzianoski, M. F. Juette, G. L. Beane, and J. Bewersdorf, “Experimental characterization of 3D localization techniques for particle-tracking and super-resolution microscopy,” Opt. Express 17, 8264–8277 (2009). [CrossRef] [PubMed]

Accuracy of the Gaussian Point Spread Function model in 2D localization microscopy

Abstract

1. Introduction

2. Vectorial theory of imaging of a dipole emitter

2.1. Preliminaries and definitions

2.2. Model for the PSF for arbitrary aberrations

2.3. Orientational averaging

2.4. The case with azimuthal symmetry

3. The effect of aberrations and dipole orientation on spot shape and localization accuracy

4. Numerical results

4.1. Effect of dipole orientation

4.2. Effect of aberrations

4.3. Effect of background

5. Conclusion

A. Temporal and ensemble orientational averaging

References and links

Cited By

Figures (6)

Equations (39)

Optics Express