We demonstrate three-dimensional tracking of fluorescent microparticles, with a computational optical system whose point spread function (PSF) has been engineered to have two twisting lobes along the optical axis, generating a three-dimensional (3D) double-helix (DH) PSF. An information theoretical comparison in photon limited systems shows that the DH-PSF delivers higher Fisher information for 3D localization than the standard PSF. Hence, DH-PSF systems provide better position estimation accuracy. Experiments demonstrate average position estimation accuracies under 14nm and 37nm in the transverse and axial dimensions respectively. The system determines the 3D position of multiple particles with a single image and tracks them over time while providing their velocities.
©2008 Optical Society of America
Discrete fluorescent particles are often encountered in applications like biological imaging [1, 2], single molecule imaging [3, 4], and flow cytometry . In most cases, these fluorescent particles move in all three spatial dimensions with time. Accurately detecting these position variations is critical in many applications. Conventional fluorescent microscopes are good for detecting position variations of focused particles along the transverse dimensions . However, because of the finite depth of field associated with the point spread functions (PSFs) of these microscopes, small axial position variations of focused particles escape undetected. Confocal microscopes offer better axial sensitivity, but the motion of the fluorescent particles precludes the use of mechanical scanning to obtain three-dimensional (3D) position information. 3D tracking requires that the 3D positions of particles at a particular time be instantaneously determined, for example by recording a single image at a given time.
One way to increase axial position estimation sensitivity is to slightly defocus the system to operate just outside focus . This increased axial sensitivity, however, comes at a price. Because of the lowering of signal levels, defocused systems have reduced transverse position sensitivity than focused systems. In photon-limited applications, defocusing buries the signal in noise, and is consequently not the best approach for 3D tracking. A weak cylindrical lens in a microscope’s imaging path creates an astigmatic PSF, which is often used for 3D tracking [8–10]. Unfortunately, astigmatic PSFs expand in size with defocus, resulting in reduced signal levels. Another interesting way to achieve 3D imaging without scanning is by using digital holography [11,12]. In general, digital holography suffers from inter-particle cross talk, coherent noise effects, and requires high space-bandwidth devices while demanding significant post processing. Although originally developed for coherent systems, it has recently been extended to incoherent emitters like fluorescent particles . However, this incoherent extension requires multiple images for determining the 3D particle locations, and is hence not suitable for 3D tracking.
In another report , we demonstrated single-image 3D localization of discrete scattering particles in a photon-unlimited monochromatic wide-field microscope presenting a double-helix (DH) PSF, which exhibits two lobes that rotate continuously with defocus. An information theoretical analysis showed that the DH-PSF is fundamentally better suited than the standard diffraction-limited PSF for 3D position localization in photon-unlimited systems with either Gaussian or Poisson noise.
In this article, we report our studies on photon-limited DH-PSF systems. In order to compare the 3D position localization accuracies of systems with different PSFs, we introduce a new quality measure based on the average Fisher information over the 3D volume of interest. Using this measure, we show that in the photon-limited case, DH-PSF systems carry higher average Fisher information than the standard PSF systems. We then demonstrate simultaneous 3D position estimations of multiple fluorescent particles using a single image with accuracies in the nanometer regime. By acquiring multiple images at periodic intervals, we track moving fluorescent particles in three dimensions, and calculate their velocities.
1.1. Double-helix point spread function (DH-PSF)
A DH-PSF system can be implemented by introducing a phase mask in the Fourier plane of an otherwise standard imaging system. The phase mask is designed such that its transmittance function generates a rotating pattern in the focal region of a Fourier transform lens [15–18]. Specifically, the DH-PSF exhibits two lobes that spin around the optical axis as shown in Fig. 1(a). Note that DH-PSF displays a significant change of orientation with defocus over an extended depth. In contrast, the standard PSF presents a slowly changing and expanding symmetrical pattern throughout the same region [Fig. 1(b)].
While analytical solutions for helical beams provide valuable insight on wave propagation  and can be used in photon-unlimited applications , they do not provide the high-efficiency transfer functions required for photon-limited systems. Hence, we use a design that confines the helical pattern to a specific axial range of interest to attain high efficiency systems . Unlike standard and astigmatic PSFs, the DH-PSF concentrates its energy in its two main lobes throughout this range of operation, and is consequently well suited for photon-limited applications.
2. Cramer-Rao bounds in photon limited systems
The position estimation accuracy of a PSF in the presence of noise can be quantified by computing its Cramer-Rao bound (CRB). The CRB of a PSF represents the lowest possible position estimation variance that can be achieved by an unbiased estimator based on that PSF [14,16,19–21]. CRBs are computed for different noise conditions by appropriately choosing the noise distribution and the noise level relative to the signal level. For the 3D position estimation problem, CRBs for the X, Y, and Z dimensions are obtained from the diagonal elements of the inverse of the 3x3 Fisher information matrix I [14, 19], which is calculated as follows:
where θ=[X, Y, Z], the indices m and n are either 1, 2, or 3, E is the expectation, and pi,j(k|θ) is the probability density function for the pixel in i th row and j th column.
For a given noise level, the position localization accuracy of a PSF is best when the intensity of the PSF spans the dynamic range of the detector. Because the energy of the rotating PSF is distributed in two of its main lobes, for a given photon budget and NA, the peak intensity of the in-focus rotating PSF is lower than that of the standard PSF. In order to compare the best possible accuracies of the rotating and standard PSFs, the peak intensities of both PSFs should span the detector’s dynamic range. This is possible in photon-unlimited systems, such as bright-field microscopes, where increasing the illumination intensity increases the intensity of the detected PSF.
However, in photon-limited systems, the detected PSF intensity cannot be arbitrarily increased. Here we call photon-limited systems those in which the illumination intensity cannot be arbitrarily increased and those in which even an increase in illumination intensity will not increase the detected intensity.
Tracking of fluorescent particles is indeed a good example of such a photon-limited problem. Because the excited state lifetime of a fluorescent molecule limits its emission intensity, any increase in the excitation intensity beyond a certain point would not help increase the number of emitted photons per unit time. Increasing the detector’s exposure time is a common alternative. However, in fluorescent particle tracking situations where the particles are moving, large exposure times produce motion blurs in the detected image. Such systems, therefore, have an inherent limited photon budget for every detected image.
Introducing a DH-PSF mask in a photon-limited system distributes the available photons into two main lobes. If the detector parameters (exposure time, gain) are set so that a standard PSF spans the detector’s dynamic range, the DH-PSF will not span the detector’s dynamic range because of its lower peak intensity. This fact is critical in the CRB calculations and fundamentally different than in the photon-unlimited case .
Many well-designed systems can be considered shift-invariant, i.e. the PSF does not change in shape when an object point is moved in the transverse (X, Y) dimensions. Consequently, the PSF CRBs do not vary for different transverse locations. However, the system’s PSF does change if the object point is displaced axially. Therefore, the PSF CRBs do vary for different axial locations. For a given noise level, the CRBs for X, Y, and Z dimensions for the entire 3D object volume can hence be completely described by three 1D plots: CRB X, CRBY, and CRBZ, as functions of axial distance.
Accordingly, Fig. 2 compares the DH-PSF CRB and the standard PSF CRB for a photon-limited system with Poisson noise. This model can be used for example in a photon-limited system with a shot-noise limited detector. In Fig. 2, the signal to noise ratio (SNR) of the standard PSF is 30. SNR is defined as the ratio of the in-focus peak intensity to the noise standard deviation. For a given point source power, the peak of the in-focus DH-PSF is 7.18 times lower than the peak of the standard PSF, so its SNR is also lower. However, in the Poisson case, the noise standard deviation decreases with a decrease in the signal. Hence, the Poisson DH-PSF SNR is .
For the Z dimension [Fig. 2(a)], the DH-PSF has lower CRB than the standard PSF, except in a region just outside the depth of field, where the standard PSF CRB is slightly lower. CRB values are very high in the depth of field of the standard PSF, indicating that axial position estimations in its focal region will have very poor accuracies. For the X dimension [Fig. 2(b)], the DH-PSF CRB is lower than the standard PSF CRB except exactly in the focal region. The standard PSF has identical CRBs for the X and Y dimensions because of its transverse symmetry but the DH-PSF CRB is higher for the Y than for the X dimension. The standard PSF CRB for Y dimension is lower than that of the DH-PSF in the focal region [Fig. 2(c)]. However, outside the depth of field region, the DH-PSF has lower CRB.
Since most practical systems are required to locate and track particles over an extended axial range, it is useful to quantify the information content of the PSFs in all three dimensions and capture the essence of each PSF design in a single parameter. Therefore, we introduce a metric defined as the sum of the CRBs in all three dimensions, CRB3D=CRBX+CRBY+CRBZ, and its average over the axial region of interest,
to quantify the overall performance of PSFs throughout the 3D volume.
For the parameters used in Fig. 2, the CRB3D largely mimics CRBZ, because for a 0.45NA objective, CRBZ is much higher than CRBX and CRBY. The DH-PSF performs better in the regions A and C of the plot, while the standard PSF does better in region B. CRBAVG of the standard PSF is infinity, and that of the DH-PSF is 239 nm2, which corresponds to an average combined standard deviation of 15.5nm in all three dimensions for a 0.45NA objective. DH-PSF is consequently better suited than the standard PSF for photon limited systems with Poisson noise. Moreover, the CRB3D for the DH-PSF is more uniform throughout the whole range implying that a high 3D accuracy can be achieved over a long axial range.
2.1. Practical DH-PSF estimator
While the CRB gives the lowest possible estimator variance, it does not provide the means to achieve this variance. Maximum likelihood estimators can reach the CRB but are computationally intensive because they require iterative operations. The DH-PSF estimator used here is non-iterative as it directly uses the PSF lobes for 3D position estimation. Specifically, the transverse (X, Y) locations are estimated from the midpoint between the centroids of the two PSF lobes, while the axial location is estimated by mapping the angular orientation of the line joining the lobe centroids to axial distance using a look-up table. This estimator is fast when compared to a maximum likelihood estimator but it does not reach the CRB.
We now describe a DH-PSF system (Fig. 3) designed for fluorescent particle tracking in three dimensions. A phase-only spatial light modulator (SLM) was used to encode the DH-PSF mask.
The sample consisted of 1µm wide yellow-green fluorescent microspheres with 505nm excitation and 515nm emission peaks. These microspheres were excited with the 488nm Argon laser line. A half wave plate and a polarizer produced horizontally polarized light, as required for the SLM. A rotating diffuser destroyed the spatial coherence of the beam. While this was not critical for fluorescence imaging, it was required to avoid coherent artifacts when operating the system in bright field mode. Two lenses, one with 85mm focal length and the other with 100mm focal length, collected the scattered light from the diffuser and focused it on the sample mounted on a piezo stage.
The imaging path can be separated into a microscope section and a signal processing section. In the microscope section, a 1.3NA infinity corrected oil-immersion objective and a tube lens with 140mm focal length produced an 85x magnified image of the sample. A 515nm interference filter with a 10nm bandwidth, located in the infinity space of the microscope section, blocked the excitation beam. In the signal processing section of the imaging path, two achromats with focal length 125mm were aligned in a reflective 4f configuration such that the SLM was located at a distance of 125mm from both of the achromats. The plane of the SLM was slightly tilted relative to the optical axis to avoid the need for an imaging path beam splitter, which would waste 75% of the emitted photons. While this tilt causes the SLM plane to be tilted with respect to the Fourier plane of the signal processing section, the effect of this is negligible and can be corrected. Finally, an Andor iXon detector, placed 125mm from the second achromat, recorded the image. A linear phase was added to the DH-PSF phase mask, and the combination was phase wrapped (see top-right inset in Fig. 3), before loading the mask into the SLM. This avoided the non-ideal 0th order response of the SLM by pushing the DH-PSF image to the SLM’s 1st order. The undiffracted 0th order was essentially the same as a standard PSF image as it was unaffected by the DH-PSF mask.
In order to characterize the system’s position localization accuracies, we imaged four fixed fluorescent microspheres embedded inside the volume of cured optical cement that has a refractive index 1.507. While the standard PSF image blurs out the microspheres outside the focal region, the DH-PSF image encoded the axial position of the microspheres in the PSFs’ rotating angles (bottom insets in Fig. 3). Each microsphere’s transverse and axial positions were estimated using the practical DH-PSF estimator (Section 2.1). The appendix presents details of the calibration and correction methods.
We estimated the position of the four microspheres by recording 100 successive images. From these data we were able to determine the standard deviation of the position estimation. The single-image standard deviations were (σX,σY,σZ)=(14nm,13nm,37nm), on an average among the studied particles (see appendix for details). Because these microspheres were fixed, we could average the 100 estimations. The standard deviation of the average estimate was (σx̄,σȲ,σZ̄)=(3nm, 3nm, 6nm) on an average for the four particles.
It is important to note that these accuracies do not represent the fundamental limit of the DH-PSF system in Fig. 3. CRB calculations of a similar system with 1.3NA, 100x magnification, 4.4µm pixel width, 514.5nm emission wavelength peak, and an SNR of 11.2 (as in the abovementioned example) revealed the single-image position localization limits to be (σX,σY,σZ)=(0.5nm,0.9nm,0.7nm) for an in-focus particle. These accuracies can potentially be achieved using more complex estimators, carefully calibrated detectors, and a more stable optical setup1.
3.1 Tracking experiment
We extended the procedure to track moving fluorescent particles by recording DH-PSF images at regular time intervals.
Figure 4 shows movies of 3D tracking of four fluorescent microspheres with the DH-PSF system. These microspheres are moving inside a blob of cured optical cement that has water trapped within it. While microspheres 2, 3, and 4 are moving in water, the microsphere 1 is relatively immobile as it is attached to optical cement. Except for particle 4, which moves out of the field of view and later reappears, all other particles are within the field of view for the entire time of the movie (10 seconds). 100 images of the moving microspheres are recorded with an exposure time of 100ms.
In the standard PSF movie [Fig. 4(a)], the transverse motion of the microspheres is detectable because of the corresponding transverse displacement of the PSF. However, the axial displacements of the particles cause their PSFs to blur, and consequently bury in noise. In contrast, in the DH-PSF movie [Fig. 4(b)], both transverse and axial displacements of the microspheres are detectable from their corresponding PSFs. When a micosphere moves in the transverse dimension, as in the standard PSF case, the DH-PSF cross-section also moves in the transverse dimension. Consequently, the instantaneous transverse position of the microsphere can be estimated from the midpoint of its PSF lobes. As the microsphere moves in the axial dimension, its PSF rotates. The microsphere’s instantaneous axial position can thus be determined by mapping the rotation angle of its PSF to the axial dimension using the calibration plot [Fig. 6(a)]. With the knowledge of each microsphere’s 3D locations every 100ms, the motion of the microspheres can be tracked, as shown in the Fig. 4(c) movie. Figs. 4(d-f) movie shows the X-Y, X-Z, and Y-Z projections of the 3D tracking movie of Fig. 4(c). X-Y [Fig. 4(d)] and Y-Z [Fig. 4(f)] tracking movies show that the particle 1 is immobile in the transverse dimensions, but does move a bit in the axial dimension. This subtle axial movement is undetectable in the standard PSF movie of Fig. 4(a).
3D tracking also facilitates the determination of velocities of the fluorescent microspheres. Fig. 5 shows the velocities averaged over 100ms for the particles 1 [Fig. 5(a)], 2 [Fig. 5(b)], and 3 [Fig. 5(c)], in all three dimensions.
We have demonstrated photon-limited tracking of fluorescent particles using a computational optical system with a DH-PSF. Because it is possible to estimate the 3D position of multiple particles over a long axial range with a single DH-PSF image and using a fast estimator, this technique is attractive to track multiple particles in real time in three dimensions. An information theoretical analysis shows that the DH-PSF carries higher average 3D information over a long range than a standard PSF. It also shows that subnanometer accuracies are achievable with the DH-PSF using improved estimators.
Calibration and detailed estimations
A calibration plot that maps rotation angles to axial locations was obtained by translating the piezo-stage holding the sample in steps of 0.8µm over a 5.6µm range, and by estimating the DH-PSF angle of particle 4 (in bottom inset of Fig. 3) from 20 successively recorded images at each location [Fig. 6(a)].
Because of system aberrations, including sample induced aberrations, the transverse location of an object point is not exactly the midpoint of the centroids of the two PSF lobes for all axial locations. Corrections were applied to the midpoint to determine the exact transverse location. The X and Y coordinates of the midpoint of the lobe centroids of particle 4 were plotted as a function of the axial distance to obtain correction factors, Xcr [Fig. 6(b)] and Ycr [Fig. 6(c)]. These correction factors were subtracted from the lobe midpoint location to obtain the transverse position of a particle. Similar curves for other particles within the field of view confirmed that these corrections are essentially shift-invariant. These plots also account for an additional minor correction that compensates for the slight asymmetry of the DH-PSF (see Fig. 1), which also shifts the midpoint of the lobes’ centroids as the system is defocused. For the system in Fig. 3, the standard deviation of this movement was theoretically calculated to be 2nm and 5nm for the X and Y dimensions.
Table 1 lists the positions (X̄,Ȳ,Z̄) and standard deviations for all four particles in Fig. 3. σX,σY, and σZ refer to the standard deviations along the X, Y, and Z dimensions, respectively. σX and σY of a particle were computed from the standard deviation (σXm,σYm) of the midpoint (Xm,Ym) of the lobes of the particle, and the standard deviation (,) of the correction factors averaged over the 20 measurements. σZ is computed from the standard deviation (σZa) of the rotation angle (Za) estimation, and the standard deviation () of the Z calibration (Zcl) averaged over the 20 measurements. Specifically, , , and . σXm,σYm, and σZa are obtained as the standard deviation of the 100 estimates. , and are calculated as the standard deviation of the average of 20 estimates.
The standard deviation of the average of N estimates is smaller than the standard deviation of each of N estimates by a factor of (N)1/2. Because the correction and calibration plots are obtained from an average of 20 estimates, the standard deviations ()of the average estimates are determined by dividing the standard deviation of each estimate by (20)1/2. However, σXm, σYm, and σZa are the standard deviations of a single estimate. σX, σY, and σZ listed in the table therefore represent single-image position estimation standard deviations. In other words, σX, σY, and σZ are the accuracies when only one image is used for position estimation.
For fixed particles, the estimation accuracy can be improved by recording multiple images. Indeed, when estimates of Xm,Ym, and Za from 100 images are averaged, the average estimates (, and ) will have a factor of (100)1/2 reduction in standard deviation (). The standard deviations of the average position estimates (X̄,Ȳ,Z̄) are , , and . The average values of σX̄,σȲ, and σZ̄ for all four particles were 3nm, 3nm, and 6nm, respectively.
It is worth emphasizing that (typically slower) optimal estimators and tight vibration control conditions could help reach the DH-PSF CRB limit for a given SNR.
|1||In addition to the photon noise considered in the CRB calculations, the above experimental standard deviations are also affected by the nanoscale vibrations in the optical setup caused by components such as the water-cooled Argon laser and the rotating diffuser. Tighter vibration control conditions could help approach the DH-PSF CRB limit for a given SNR.|
S. R. P. Pavani thankfully acknowledges support from a CDM Optics - OmniVision Technologies Ph.D. fellowship. This work was partially funded by the Technology Transfer Office of the University of Colorado and the National Science foundation (award ECS-225533 and DMI-0304650).
References and links
4. G. J. Schütz, M. Axmann, and H. Schindler, “Imaging single molecules in three dimensions,” Sing. Mol. 2, 69–74 (2001). [CrossRef]
7. M. Speidel, A. Jonas, and E.-L. Florin, “Three-dimensional tracking of fluorescent nanoparticles with subnanometer precision by use of off-focus imaging,” Opt. Lett. 28, 69 (2003). [CrossRef] [PubMed]
8. H. P. Kao and A. S. Verkman, “Tracking of single fluorescent particles in three dimensions: use of cylindrical optics to encode particle position,” Biophys. J. 67, 1291–1300 (1994). [CrossRef] [PubMed]
9. L. Holtzer, T. Meckel, and T. Schmidt, “Nanometric three-dimensional tracking of individual quantum dots in cells,” Appl. Phys. Lett. 90, 053902 (2007). [CrossRef]
11. S. -H. Lee, Y. Roichman, G. -R. Yi, S. -H. Kim, S. -M. Yang, A. van Blaaderen, P. van Oostrum, and D. G. Grier, “Characterizing and tracking single colloidal particles with video holographic microscopy,” Opt. Express 15, 18275 (2007). [CrossRef] [PubMed]
13. J. Rosen and G. Brooker, “Non-scanning motionless fluorescence three-dimensional holographic microscopy,” Nature Photon. 2, 190–195 (2008). [CrossRef]
14. S. R. P. Pavani, A. Greengard, and R. Piestun, “Three-dimensional localization with nanometer accuracy using a double-helix point spread function system,” (submitted).
17. Y. Y. Schechner, R. Piestun, and J. Shamir, “Wave propagation with rotating intensity distributions,” Phys. Rev. E 54, R50–R53 (1996). [CrossRef]
18. R. Piestun, Y. Y. Schechner, and J. Shamir, “Propagation-invariant wave fields with finite energy,” J. Opt. Soc. Am. A 17, 294–303 (2000). [CrossRef]
19. S. M. Kay, Fundamentals of Statistical Signal Processing: Estimation Theory (Prentice-Hall, 1993).
20. T. M. Cover and J. A. Thomas, Elements of Information Theory (Wiley-Interscience, 1991). [CrossRef]