2.1.

Spiral Phase Plate

A perspective view of the SPP is shown in green in Fig. 1. The central SPP design idea is to generate a coherent superposition of several nondiffracting Bessel vortex modes using a SPP placed in the aperture of an imaging system. Compared with LG modes,^30–32 the Bessel modes^33,34 have significantly reduced diffraction as they propagate. In our case, the SPP employs seven concentric annular zones, or rings, of spiral phase retardation with changing orbital angular momentum (OAM) quantum number between the zones. Designing the SPP so that this OAM quantum number increments in regular steps between successive zones and the initial phase discontinuity lines of all zones are azimuthally aligned allows the Bessel beams to interfere constructively and yield a compact image that, for a specific choice of zone radii, rotates uniformly and unidirectionally with changing depth of the light emitter. An approximate analytical expression for the ROR of the emitter image with its depth change may be derived from Eqs. (18) and (19) of Ref. 9 using the latter equation to express the left-hand side of the former in terms of the ratio, $l/L$, and then setting it equal to the spiral phase, $l\theta$, corresponding to the $l$’th zone. Differentiating that equality yields the ROR, which may be expressed in terms of the system parameters as

Fig. 1

Schematic diagram for encoding 3D source location in 2D images. Two point sources (orange and green) are shown at the left. The orange source is in axial focus. The green source is displaced along $x$ and $z$ with $\delta z<0$. The seven zones of the SPP each impart a unique OAM to the imaging wavefront with the OAM eigenvalue of the $l$’th zone being $l\hslash$. A specific choice of the zone radii allows the seven OAM wavefronts, two of which are shown spiraling, to interfere constructively at the camera plane to create a bright spot that is rotated about the Gaussian image of the source by an angle proportional to the source misfocus ($\delta z$). See text for details.

The original SPP design described by Prasad⁸ has a number of concentric annular zones, each with a continuous spiral phase ramp. The $l$’th zone consists of a single azimuthal ramp (vortex) such that the optical phase retardation increases continuously from zero to $2\pi l$ radians, whereas the azimuthal angle goes from zero to $2\pi$ radians, requiring an optical phase retardation as large as $14\pi$ radians, for example, for a seven-zone design. Manufacturing such an SPP to the required precision is not possible with current DOE fabrication technology. This major practical obstacle inherent in the original SPP design was circumvented by Ritsch-Marte’s group¹¹ in 2014 by subdividing the $l$’th annular zone into $l$ subramps, each of which increases in optical phase retardation from 0 to $2\pi$ at the design imaging wavelength, like a spiral staircase (Fig. 1). The present paper is a study of the performance of a fused silica DOE using such a modified phase-wrapped design with further modifications included to account for the high NA of modern microscopes.⁹ Our high-NA SPP was fabricated by Vortex Photonics, Planegg, Germany.

2.2.

Microscope Setup and Image Analysis

An Olympus IX-71 microscope was used for the experiments with fluorescent microbeads. The objective had a magnification of $60\times$, NA = 1.3 and used silicone oil (refractive index = 1.406) as the immersion medium. Images were taken under two scenarios: a fixed sample stage and a moving objective or a moving sample stage and a fixed objective, as shown in Fig. 2(a). In the first scenario (moving objective), movement along the $z$-axis was controlled using a Physik Instrumente PIFOC objective scanner (P/N PI72Z2CAQ); in the second scenario (moving stage), the $z$-position of the sample was controlled by a piezoelectric Mad City Labs stage (Model MCL-MOTNZ). In both cases, the relative position between the sample bead plane and the objective was moved in 100-nm steps from $-3$ to $+3\mu \mathrm{m}$ relative to the sample focal plane. This corresponds to 61 positions for each of which an image was recorded. Images were taken without and with the SPP in place. The SPP was in a standard nosepiece holder from Olympus that fit directly under the exit of the objective, which could easily be slid into or out of place.

Fig. 2

(a) Schematic diagram of the microscope setup, shown for an upright system. (b) SPP structure viewed from above, and falsely colored to show height of structures. Colormap scale is in μm with the height range from 0 to 1.13 μm. (c) Glancing view of actual SPP showing the structure of the surface. (d) SPP mounted in the Nomarski prism slider, for insertion in the microscope nosepiece. (b) and (c) Courtesy of Vortex Photonics, Planegg, Germany.

Bead samples were illuminated by an LED with a peak wavelength of 470 nm (Thorlabs model 470L4). The detector was a Hamamatsu ORCA Fusion sCMOS camera with exposure set to 200 ms. Standard procedure was to select a region of interest (ROI) of $500\times 500\text{pixels}$ ($54\times 54\mu {\mathrm{m}}^2$), containing 4 to 20 beads. The LED intensity was adjusted to fill the intensity range of the 16-bit camera. The maximum intensity occurred in the middle plane (focal plane depth within the slide) where the hotspot of the RPSF pattern was most compact and most intense. The fluorescent beads had SNRs at the focal plane of $\sim 100$ and at extreme depth detunings of $\pm 3\mu \mathrm{m}$, the SNR was 28. Both SNR values were calculated using the maximum intensity of the hotspot compared with nearby background levels. The SNR values for cell samples, on the other hand, were 27 and 19, respectively, for the first and last frames in the time series, where photobleaching is the cause of the decreasing fluorescent intensity with time.

Ten replicates of the same field of view (FOV) were collected to improve the SNR of the images. These replicates were automatically collected using micromanager and Pycro-Manager software to control the motion along $z$ of the stage as well as the camera image recording. A total of 11 FOV’s representing 153 distinct beads were analyzed.

A similar setup, based on a motorized Olympus IX-83, was used for live cell imaging. Position in $x/y$ was controlled using an ultrasonic motorized stage (Olympus BX3-SSU), whereas the $z$ position was controlled using the motorized nosepiece, enabling collection of $z$-stacks with 100-nm steps over a $6\text{-}\mu \mathrm{m}$ range, as for the bead samples. Sample drift was avoided using the $z$-drift compensation system (Olympus IX3-ZDC). Cells were kept in physiological conditions with a stage-top incubator (Tokai Hit), maintaining 37°C and 5% ${\mathrm{CO}}_2$. Images were recorded with a Hamamatsu ORCA Flash 4.0 camera, controlled by the CellSens software (Olympus).

2.3.

Sample with Fluorescent Beads

We used an Olympus objective designed specifically to match the optical environment of a cell. Thus, our fluorescent bead sample was optically designed to match fluorescent emitters inside a cell. To achieve this, the sample consisted of a microscope slide coated with a flat plane of 100-nm-diameter yellow-green fluorescent polystyrene microspheres (Fluoresbrite YG microspheres, Polysciences Inc, Warrington, Pennsylvania, United States) embedded in a plane between two layers of optical adhesive (Norland Products #NOA 139, Jamesburg, New Jersey, United States), whose refractive index, 1.39, closely matched the optical environment of cells.

2.4.

Mammalian Cells

U2OS osteosarcoma cells with a stable genomic integration of Lactose repressor (LacR) binding arrays (256 copies of Lactose operon - LacO³⁵) were used to track chromatin motion in 3D. The cells were seeded in 35-mm glass-bottom dishes (MatTek) at 100,000 cells per dish and maintained in Dulbecco’s modified Eagle medium (DMEM) supplemented with 10% FBS (Sigma) in a humidified incubator set to 37°C, 5% ${\mathrm{CO}}_2$. To visualize the LacO repeats, LacR fused to green fluorescent protein - GFP (LacR-GFP) was transiently expressed by transfection of the corresponding plasmid DNA using Lipofectamine 3000 (Invitrogen), 24 hr before imaging. For fixation, the cells were incubated for 20 min in 4% paraformaldehyde solution (Sigma cat# HT5011).

2.5.

Template Matching to Determine 3D Position from 2D Images

To recover the $x$, $y$, and $z$ positions of the beads from the 2D camera images, we used a library of 61 images (templates), each being the image of a single bead. The images were acquired with the SPP in place. The $z$ position of the beads (defocus) fell between $-3$ and $+3\mu \mathrm{m}$ about the plane of Gaussian focus. The template images were $51\times 51\text{pixels}$ with the bead centered in each image. For each of the 10 image sets in a FOV, a template bead was selected to create the template library. The template bead was the bead that produced the most symmetric pattern for the average intensity $z$-projection of all 61 planes. A study of many beads in multiple FOVs indicated that this optimal point of minimal spherical aberration was located about $5\mu \mathrm{m}$ to the left and $10\mu \mathrm{m}$ above the actual center pixel of the camera FOV. Planes were collected so that plane #1 was at or close to $\delta z=-3\mu \mathrm{m}$, and plane #61 was close to $\delta z=+3\mu \mathrm{m}$.

Template matching was applied in two stages. In the first stage, we used a built-in MATLAB function, normxcorr2, that computed the intensity-based normalized 2D cross-correlations between each template image and each target image. This function translates the “moving” target image one pixel at a time in both $x$ and $y$ directions. The criterion for best match is the maximum in the 2D cross-correlation value. This method is limited to pixel precision, $\pm 0.5\text{pixel}$, or about 54 nm, in $x$ and $y$. In the second stage, we used paraboloid surface fit technique³⁶ to find $x$ and $y$ to subpixel precision. For the determination of $z$, the second stage consisted of using a simple one-dimensional (1D) quadratic curve fit to the three points around the optimal template for a given target plane. In this way we could find the $z$-axis localization to subpixel precision.

2.6.

Tracking Fluorescently Labeled DNA Regions in Live Cells

The template-matching approach (Sec. 2.5) was used on live cells as in the case of fluorescent beads. A template library of images was generated from 61 planes taken for a fluorescently labeled Lac array in a fixed cell. The Olympus IX-83 motorized stage was moved in 100-nm increments over a $6\text{-}\mu \mathrm{m}$ range while collecting fluorescent images. A $51\times 51\text{pixel}$² ROI was created around the center of the fluorescent region in the cell nucleus. We also used the MOSAIC Particle Tracker (ImageJ plugin)³⁷ to track the same fluorescently labeled DNA spots with the SPP removed, i.e., with no SPP. This allowed us to compare our 3D tracking results with classic 2D tracking.

2.7.

Determining the Diffusion Coefficient from Tracking Data in Live Cells

Mean-squared displacement (MSD) values were computed for all time points ($t$) according to the following equation:

2.8.

Vector-Field Simulation of the RPSF for a Microscope with a Single-Spot Generating SPP

Our vector field simulations of the rotating PSF were based on expressions in Ref. 9. The appropriate equations were translated into MATLAB code, which was used to generate the theoretical RPSFs shown in Fig. 3(a) and to make predictions about our experimental observations. The MATLAB code is available as a Supplementary Material in Code 1 (Ref. 38).

Fig. 3

(a) Montage of theoretical spiral patterns, represented in false color (see the color bar for values of the image intensity), produced by the single-spot SPP both above and below the nominal focal plane for the same optical parameters as those used to collect experimental data. The value of the defocus $\delta z_0$, in μm, is indicated in parentheses in each image. (b) Montage of the corresponding experimental patterns. The maximum intensity in each sub-image in both montages is normalized to 1. Panels (a) and (b) are represented at the same scale, as they would appear projected into the sample.

3.1.

Optical Performance of the SPP

We recorded 61 RPSFs of a 100-nm-diameter spherical fluorescent bead with the source bead located in the range $\delta z_0=-3$ to $+3\mu \mathrm{m}$. The Gaussian focal plane is defined by $\delta z_0=0\mu \mathrm{m}$. Figure 3(a) shows a montage of seven theoretical images of a 100-nm bead with $z$-gaps of $1\mu \mathrm{m}$ (total range $-3$ to $+3\mu \mathrm{m}$) generated by vector field calculations using the design parameters of the SPP with seven zones, as shown in Code 1 (Ref. 38). The beads were endowed with a Lambertian reflectance profile in our simulations. The length scale matches the scale expected under the experimental conditions of the microscope system. Intensities are normalized to 1 for each subimage. Normalized experimental sub-images collected for the same seven planes are shown in Fig. 3(b). Similarities between the two sets are the rotations of the maximum intensity regions with depth, the creation of spiral arms of growing width as depth detuning increases, and the reversal of helicity of the spiral pattern as one goes from negative detunings ($z<0\mu \mathrm{m}$) to positive detunings ($z>0\mu \mathrm{m}$) and vice versa. There are also qualitative differences between the calculated and observed “spiral galaxies” outside the “hot spot,” including a broader, more diffuse footprint of the experimental images at larger distances from focus.

There are several key quantitative metrics to characterize the optical behavior of the SPP. The original goal motivating the design was to produce a single spot rotating on a circle with the angle of rotation depending linearly on depth. The selection of the zones (i.e., the number, type, and area of each zone of the SPP) was based on these criteria.

As the plane is moved along $z$, the hotspot from a bead rotates in a circle about the $z$ axis (Figs. 1 and 3). The bright spot ideally moves along a circle of fixed radius as the detuning from focus changes. We identified the centers of the bright spots at each position, using a center of mass intensity algorithm. Figure 4(a) shows the trajectory of the intensity peaks in the theoretical images [like those shown in Fig. 3(a)] for a range of depths corresponding to $-2$ to $+3\mu \mathrm{m}$ in increments of 100 nm. The range from $-3$ to $-2\mu \mathrm{m}$ is not shown because most of the experimental datasets in this range had their spots split into two angularly separated spots, preventing us from determining spot center positions. With the design of the SPP, the theoretical points from $-2$ to $+3\mu \mathrm{m}$ fill an angular range that is about 1.2 radians ($\sim 70\mathrm{deg}$) $>2\pi$. In Fig. 4(b), we show corresponding experimental data with the same range of depths and increments. The 51 measurements follow the theory well. The circle fits to these trajectory plots allow us to compare the experimental results to theory. The predicted radius of the trajectory circle in the camera plane was $39.7\mu \mathrm{m}$. The observed radius for the specific bead in Fig 4(b) was $39.9\mu \mathrm{m}$. The measured radius from the combined data on 56 beads was $39.1\pm 1.4\mu \mathrm{m}$ [mean ± standard deviation (SD)]. The theoretical and experimental values for the circle radius differ by 2%.

Fig. 4

Hotspot trajectories: theory versus experiment. (a) Theoretically calculated hotspots from a set of 51 images of the irradiance pattern corresponding to depths in the range $-2$ to $+3\mu \mathrm{m}$ in intervals of $\delta z=100\mathrm{nm}$. (b) Corresponding experimental measurements. Also shown (in blue) is the best fit circle to the hotspot data. The first point at $z=-2\mu \mathrm{m}$ is indicated by a black dot. Scale is that in the camera plane.

Another important quantitative metric is the design criterion of having a uniform ROR, i.e., the hotspot should rotate on a circle such that the angle versus $z$-depth plot is linear. An approximate analytical expression for the rotation rate of the bright spot in the intensity distribution produced by a point emitter is given in Eq. (1). This analytical rotation rate depends on several physical parameters, such as the objective $\mathrm{NA}=1.3$, the center wavelength of the emitter, ${\lambda }_0=525\mathrm{nm}$, the largest OAM quantum imparted to photons by the SPP ($L=7$), and the refractive index of the medium, $n=1.406$. Substituting these values into Eq. (1) gives the theoretically predicted result of $d\theta /dz=1.49\mathrm{rads}/\mu \mathrm{m}$. Correspondingly, an axial range of $2\pi /1.49=4.2\mu \mathrm{m}$ is indicated for one complete rotation of the pattern around the center, which underestimates the experimentally observed range, as seen in Fig. 4 and which we discuss further near the end of the present subsection.

To obtain a “theoretical” value for the expected rotation rate under high NA conditions, we processed the images under our exact high NA conditions ($\mathrm{NA}=1.3$) and then selected the coordinates corresponding to the highest intensity pixel in the computed image. We plotted the data for the 51 points in Fig. 4(a) and found the slope by performing a linear fit to the data [Fig. 5(a)]. Here the red circles correspond to the maximum values of the points in the theoretically generated images for the different depths. The blue curve is the best linear fit to the red circles and gives a resulting theoretical slope of $d\theta /dz=1.50\mathrm{rads}/\mu \mathrm{m}$, in excellent agreement with the analytical formula given in Eq. (1).

Fig. 5

Hotspot rotation versus source depth: theory versus experiment. (a) Angles versus depth along $z$ of theoretically calculated hotspots from a set of 51 images of the irradiance pattern corresponding to depths in the range $-2$ to $+3\mu \mathrm{m}$ in intervals of $\delta z=100\mathrm{nm}$. Changes of rotation angle taken from the points in Fig. 4(a) with the zero of angle arbitrarily taken to be the angle at the point $\delta z_0=-2\mu \mathrm{m}$. (b) Corresponding angles of experimentally observed hotspots (red circles) for the bead shown in Fig. 4(b). The hotspot rotates clockwise with increasing $z$ as clearly seen from Figs. 3 and 4. The blue line is the best fit line.

The corresponding experimental plot of angle versus $z$-depth is shown in Fig. 5(b) and this is taken from the trajectory plotted in Fig. 4(b). The red circles correspond to the center-of-mass of intensity values that equal or exceed 0.9 of the maximum intensity. Figure 5(b) shows the experimental data for the template bead in a FOV with four beads. Here we show the same range as in Fig. 4(b), namely $\delta z_0=-2$ to $+3\mu \mathrm{m}$. For the same set of data that we used to report the trajectory radius (trajectories of 56 beads in 56 FOVs distributed over 2 samples), we found that $⟨d\theta /dz⟩=0.99\pm 0.08\mathrm{rads}/\mu \mathrm{m}$ (mean ± SD). This experimental ROR is 34% lower than the value given by the analytical expression, Eq. (1). Note that in Fig. 5(b), the slope corresponds to that of a single bead in a single FOV, so its slope differs from the mean value above that corresponds to 56 beads in different FOVs.

One possible explanation for the discrepancy between the predicted and observed RORs is that since the SPP is not exactly at the exit pupil plane (EPP) of the objective, we would need to correct for their separation. Equivalently, one must map the SPP back to that plane. For rays diverging when incident on the SPP, which is the case in our microscope, with the SPP placed downstream in the turret nosepiece, this would mean that the SPP would appear both demagnified and diffractively blurred when mapped back to the EPP of the objective. A demagnified SPP is equivalent to a larger effective $L$ (note the zone radii are inversely proportional to the square root of $L$ for low NA with a more complicated decrease with $L$ for high NA). Changing effective $L$ from 7 to about 10 would bring the theoretically predicted ROR into agreement with the experimentally observed ROR and the theoretical maximum axial range for one complete rotation, over which the pattern deteriorates unacceptably, into agreement with the observed range of about $6\mu \mathrm{m}$.

The diffractive blurring of the SPP when mapped back to the objective EPP, on the other hand, is possibly responsible for the second important difference between the theoretically predicted and observed bead images, namely the broader, more diffuse footprint of the latter. Note, however, that despite an imperfect SS-RPSF design spurred by our desire to achieve simplicity of implementation, we still successfully tracked the diffraction-limited spots with this DOE in this configuration with good localization precision for DNA loci in live cells, as discussed in Sec. 3.3.

3.2.

SPP Localization Performance

To assess SS-RPSF localization performance, for each FOV we used a $59\times 59\text{pixel}$² ROI around the template bead. Template matching was performed in a two-stage process that first used the 2D cross-correlation function in MATLAB to determine which template image best matched each of the 61 $z$-plane ROI target images for each bead in the FOV. This template-matching method only used translation along $x$ and $y$ to find the best template to match a given bead’s ROI $z$-plane image. The minimum transverse translation step was 1 pixel.

For stage two, we found the location in the lateral dimensions to subpixel precision using a quadratic surface fitting procedure that used a $3\times 3\text{pixel}$ region around the location of maximum 2D cross correlation and fit a general paraboloid (allowing both elliptic and hyperbolic paraboloids) to estimate where the cross-correlation maximum occurred between integer pixel shifts of the image and template. A nonlinear least-squares fitting function in MATLAB was used to generate the best fitting parameters, namely the $x\text{-}y$ location that gave the maximum cross-correlation.

To locate the cross-correlation maximum along $z$ to subpixel precision, we fitted a parabola to the 3 values of $z$ below, at, and above the best fit template plane. The vertex of that parabola was determined analytically. In both the transverse and axial cases, we also computed the standard deviation of the localization differences between adjacent $z$-steps. This is a metric for uncertainty in the precision of changes in location during tracking.

Overall, this simple two-step cross-correlation template-matching approach with only translational degrees of freedom gave reasonably precise localization results, as described below in detail. Our goal was to meet or exceed the localization/tracking precision achieved in our earlier work in 2D where we localized and tracked chromatin microdomains with a precision of about 40 nm.³⁹ Here we show that the mean transverse localization precisions along $x$ and $y$, as they differ from ground truth, were $⟨{\sigma }_x⟩=21\pm 13\mathrm{nm}$ and $⟨{\sigma }_y⟩=24\pm 6\mathrm{nm}$, (mean ± SD) over the entire $6\text{-}\mu \mathrm{m}$ depth range; the mean axial localization precision difference from ground truth was $⟨{\sigma }_z⟩=70\pm 15\mathrm{nm}$ over the $6\text{-}\mu \mathrm{m}$ depth range. Considering a smaller depth range ($2.9\mu \mathrm{m}$), the localization precision values are: $⟨{\sigma }_x⟩=19\pm 10\mathrm{nm}$, $⟨{\sigma }_y⟩=24\pm 4\mathrm{nm}$, and $⟨{\sigma }_z⟩=50\pm 3\mathrm{nm}$ (Fig. 7). These results with fluorescent latex spheres match or exceed our previous 2D results for transverse localization of fluorescent chromatin microdomains in cells. Localization in the axial direction is slightly less precise, as expected. This data are based on beads that were at least $4.32\mu \mathrm{m}$ (40 camera pixels) from their nearest neighbor and included 77,410 bead images.

Our $x$, $y$, and $z$ precisions are also comparable with those achieved by other state-of-the-art 3D localization techniques such as 3D STochastic Optical Reconstruction Microscopy (STORM),¹³ double-spot PSFs,^3,11 and phase-ramp PSF,¹⁴ but our results extend over an axial depth range that is 1.5 to 3 times larger, with a simple low-cost retrofitting of a commercially available microscope. Although some more recent approaches can achieve even better $z$ precisions that are comparable to the $x$, $y$ precisions of order 20 nm, they either do so over a much smaller axial range or require highly sophisticated optical-bench setups. The use of supercritical-angle fluorescence⁴⁰ can generate such high $z$ precision but only over a small sample thickness of order 150-nm contiguous to the cover slip, while rotating PSF interferometric imaging⁴¹ can do so over 2- to $3\text{-}\mu \mathrm{m}$ depths but only by means of an expensive dual-opposed-objective layout.

Most recently, deep-learning methods, such as DeepSTORM3D,⁴² have exhibited the potential to increase the localization precision further, even with spatially highly extended PSFs, such as overlapping tetrapod PSFs.⁷ Deep-learning methods have also been successfully applied to optimize paired PSFs to achieve precise 3D localization at high densities.⁴³ Similar machine-learning approaches with the single-spot RPSF such as that employed here can also potentially improve its error performance even under lower SNR and higher emitter-density conditions.

3.2.1.

Localization agreement between ground truth and experiment (metric #1)

The key metric for evaluating how well particle motion can be tracked in a movie is the correct determination of the change in position from one frame to the next, i.e., an accurate measure of the change in position between adjacent frames.

Transverse motion

In stepping along $z$, we did not intentionally move the sample along the transverse directions. Thus, we expect adjacent $z$-planes would have the same $x$, $y$ values. We analyzed the motion to subpixel localization precision using the two-step process described above. Figures 6(a) and 6(b) show frequency (percentage) histograms of the subpixel localization precisions for motions along $x$ and $y$, respectively, compared with the ground truth of no motion (change of 0 pixels in both $x$ and $y$).

Fig. 6

Localization agreement between ground truth and experiment. (a) Histogram of the percentage of bead images whose change in $x$-axis location agrees ($\mathrm{Dx}=0\text{pixels}$ per step, blue) or disagrees ($\mathrm{Dx}\ne 0\text{pixels}$ per step) with ground truth, for axial steps. (b) Histogram of the percentage of bead images whose change in $y$-axis localization agrees ($\mathrm{Dy}=0\text{pixels}$, orange) or disagrees ($\mathrm{Dy}\ne 0\text{pixels}$) with ground truth, for axial steps in $z$. (c) Histogram of the percentage of bead images whose change in $z$-axis location agrees ($\mathrm{Dz}=1\text{pixel}$, black) or disagrees ($\mathrm{Dz}\ne 1\text{pixel}$) with ground truth, for axial steps in $z$.

Fig. 7

Localization precision (standard deviation – ${\sigma }_x$, ${\sigma }_y$, ${\sigma }_z$) for changes in $x$ (blue open squares), $y$ (orange open triangles), and $z$ (black open circles) compared with ground truth as a function of the $z$-depth.

Note that the localization data includes all 61 planes for each bead in the FOV, including the 10 of the 61 planes that do not follow the theoretical trajectory shown in Fig. 4. Yet, as these results show, such beads can still be located with reasonable transverse localization precision and so the localizable $z$-depth range is the full $6\text{-}\mu \mathrm{m}$ range over which we gathered data. This result meets our expectations qualitatively because the 10 frames at the larger depths have distinctive intensity distributions visually.

Axial motion

Because each bead was moving in 100-nm steps along $z$ between successive images (for 61 planes), we expect that the difference in the template images that matched these would differ by 1 plane or 100 nm. Thus, we performed the same two-step template matching procedure described at the beginning of this section, which gave us an estimate of the subpixel localization precision. Figure 6(c) provides the histogram frequency data for the distribution of $z$-step values measured by template-matching, where the ground truth is that $\delta z_0=100\mathrm{nm}$ for each step. One pixel in Fig. 6(c) is 100 nm. Clearly there is a peak at a change in pixel of 1 in the distribution as expected. Also, to be expected are subsidiary peaks at nearby integer values (0 and 2), but with points in between due to the use of subpixel localization precision methods (parabolic fitting in 1D). This distribution differs significantly from a Gaussian due to its high peak at 1 as well as the subsidiary peaks at 0 and 2. We note that the high peak at one is due to the naturally larger number of points that would land at the intended change, but also note that it is partially due to assigning the subpixel fraction to 0 for cases where the second step of fitting a quadratic function to the three planes (with the template-identified plane in the middle) created a result that indicated a move of one or more pixels. We assigned these cases a fractional subpixel shift of 0 since we do not allow a subpixel-determining operation to shift our final position by one pixel or more.

An estimate for the $z$-axis localization precision is given in the next section where we separately analyze the localization precision for the 60 different steps between the 61 planes through the whole range of depths.

3.3.

Tracking DNA Loci in Live Cells

We expressed the LacR-GFP and a nuclear localization signal in U2OS cells with a stable genomic integration of LacO arrays to fluorescently label subdiffraction spots of chromatin in cell nuclei [Fig. 8(a)]. This method is well established to study chromatin dynamics in yeasts and mammalian cells. LacR-GFP signals could easily be imaged using the SPP [Figs. 8(b) and 8(c)]. First, a library of cell images was generated, using fixed cells with LacR-GFP signals, moving the objective by 100-nm steps, as for the fluorescent beads sample (Fig. 9). This library was used for template-matching as described above. Note that the template library generated from fixed cell data is used for all subsequent live-cell data analyses. There was no need to generate new template libraries for each cell or biological replicate. Next, we tracked LacR-GFP spots for 60 s (300 ms exposure for 200 frames) in live cells. Figure 8(d) plots the MSD for a LacR-GFP trace, corresponding to the cell shown in Fig. 8(b). In addition to the overall MSD given by Eq. (2), we also plot the separate components ($x/y/z$) for comparison.

Fig. 8

(a) Schematic of a cell with a genomic integration of LacO DNA arrays and expressing LacR-GFP to visualize these arrays. (b) Image of a LacR-GFP fluorescent spot within the nucleus (boundary highlighted with white dotted line). Images were taken without or with SPP. (c) Enlargement of the regions delineated with yellow dotted lines in panel (b) and displayed at different time points. (d) Plot of the MSD versus time for the spot shown in panel (a). MSD for all three dimensions (blue solid circles), MSD along $z$ (green diamonds), MSD along $y$ (magenta squares), and MSD along $x$ (black open circles).

Fig. 9

Two LacR-GFP DNA spots in different U2OS cell nuclei imaged as a function of $z$-depth. The two spots are at different depths in the sample (Video 1, MPEG, 126 KB [URL: https://doi.org/10.1117/1.JBO.27.12.126501.s1]).

MSDs can be used to extract diffusion coefficients, $D$, which correspond to the fitted slope of the MSD versus $\tau$ plot divided by 6 [Eq. (2)]. Figure 8(d) shows the linear fit to the 3D motion data, which yielded a diffusion coefficient value of $D=301{\mathrm{nm}}^2/\mathrm{s}$, in agreement with previous reports on chromatin diffusion.^45,46 It also agrees with analysis of the 2D motion of the same LacR-GFP spot but tracked with the SPP removed (no SPP). In that case, we found that $D$ (no SPP) = $380{\mathrm{nm}}^2/\mathrm{s}$.

3D single-particle tracking of fluorescently labeled DNA subregions has been reported in live yeast cells by Moerner’s group using both a double-helix phase plate⁴⁷ and a tetrapod phase plate.⁴⁸ These previous studies demonstrate the important role 3D tracking can play to better understand chromatin dynamics in cell nuclei, which is the goal of our investigations. Both studies from the Moerner group involved tracking two GFP-labeled LacO/LacI foci in yeast chromosomes, near the Galactose (GAL) locus. Diffusion coefficients, derived from motions measured over 10 s (0.1-s intervals) were $\sim 2000{\mathrm{nm}}^2/\mathrm{s}$. There are differences between the time scales in this previous study (10 s) and our data (60 s). Moreover, the nuclear environments of yeast and mammalian cells differ, further contributing to different measurement outcomes. In a more recent publication, the same group studied anomalous diffusion of chromatin regions in human cells.²⁴ The DH-PSF was used to track chromatin loci in 3D, on fast (∼several minutes) and slow (hours) timescales. The effective lateral and axial chromatin diffusion values were $\sim 850$ and $\sim 700{\mathrm{nm}}^2/\mathrm{s}$, respectively, in the same range as those measured here. The differences are likely due, at least in part, to the fact that different genomic positions were tracked in the two studies.