Numerical study of the subtraction threshold for fluorescence difference microscopy

Subtraction microscopy has recently been promoted for its compactness and simplicity of enhancing spatial resolution. However, until now, the subtraction factors used in such microscopes to process raw images have been chosen experientially, and it has been impossible to determine whether the resolved structures after subtraction are appropriate or over-processed. Based on vector diffraction theory and two-dimensional convolution, this paper presents numerical investigations of the parameters that may offer the selection criterion of subtraction factors used in subtraction microscopy. It proposes two essential parameters for appropriately evaluating the subtraction factor: the fluorescence peak intensity after subtraction and resolution derivative. © 2014 Optical Society of America OCIS codes: (100.6640) Superresolution; (110.2960) Image analysis; (180.2520) Fluorescence microscopy. References and links 1. E. Betzig, G. H. Patterson, R. Sougrat, O. W. Lindwasser, S. Olenych, J. S. Bonifacino, M. W. Davidson, J. L. Schwartz, and H. F. Hess, “Imaging intracellular fluorescent proteins at nanometer resolution,” Science 313, 1642–1645 (2006). 2. M. J. Rust, M. Bates, and X. W. Zhuang, “Sub-diffraction-limit imaging by stochastic optical reconstruction microscopy (STORM),” Nat. Meth. 3, 793–796 (2006). 3. M. G. L. Gustafsson, “Surpassing the lateral resolution limit by a factor of two using structured illumination microscopy,” J. Microscopy 198, 82–87 (2000). 4. S. W. Hell and J. Wichmann, “Breaking the diffraction resolution limit by stimulated emission: stimulatedemission-depletion fluorescence microscopy,” Opt. Lett. 19(11), 780–782 (1994). 5. M. P. Gordon, T. Ha, and P. R. Selvin, “Single-molecule high-resolution imaging with photobleaching,” P. Natl. Acad. Sci. USA 101(17), 6462–6465 (2004). 6. A. Sharonov and R. M. Hochstrasser, “Wide-field sub-diffraction imaging by accumulated binding of diffusing probes,” Proc. Natl. Acad. Sci. USA 103(50), 18911–18916 (2006). 7. T. Dertinger, R. Colyer, G. Iyer, S. Weiss, and J. Enderlein, “Fast, background-free, 3D super-resolution optical fluctuation imaging (SOFI),” P. Natl. Acad. Sci. USA 106(52), 22287–22292 (2009). 8. J. Miyazaki, H. Tsurui, A. Hayashi-Takagi, H. Kasai, and T. Kobayashi, “Sub-diffraction resolution pump-probe microscopy with shot-noise limited sensitivity using laser diodes,” Opt. Express 22(8), 9024–9032 (2014). 9. J. Miyazaki, K. Kawasumi, and T. Kobayashi, “Resolution improvement in laser diode-based pump-probe microscopy with an annular pupil filter,” Opt. Lett. 39(14), 4219–4222 (2014). #223809 $15.00 USD Received 24 Sep 2014; revised 4 Nov 2014; accepted 4 Nov 2014; published 11 Nov 2014 (C) 2014 OSA 17 November 2014 | Vol. 22, No. 23 | DOI:10.1364/OE.22.028819 | OPTICS EXPRESS 28819 10. C. F. Kuang, S. Li, W. Liu, X. Hao, Z. T. Gu, Y. F. Wang, J. H. Ge, H. F. Li, and X. Liu, “Breaking the diffraction barrier using fluorescence emission difference microscopy,” Sci. Rep. 3(1441) (2013). 11. H. Dehez, M. Piche, and Y. D. Koninck, “Resolution and contrast enhancement in laser scanning microscopy using dark beam imaging,” Opt. Express 21(13), 15912–15925 (2013). 12. S. Segawa, Y. Kozawa, and S. Sato, “Resolution enhancement of confocal microscopy by subtraction method with vector beams,” Opt. Lett. 39(11), 3118–3121 (2014). 13. B. Richards and E. Wolf, “Electromagnetic diffraction in optical systems. II. structure of the image field in an aplanatic system,” Proc. R. Soc. Lond. A 253, 358–379 (1959). 14. D. Ganic, X. Gan, and M. Gu, “Focusing of doughnut laser beams by a high numerical-aperture objective in free space,” Opt. Express 11(21), 2747–2752 (2003). 15. R. Dorn, S. Quabis, and G. Leuchs, “Sharper focus for a radially polarized light beam,” Phys. Rev. Lett. 91, 233901 (2003). 16. G. M. Lerman and U. Levy, “Effect of radial polarization and apodization on spot size under tight focusing conditions,” Opt. Express 16(7), 4567–4581 (2008). 17. S. Deng, L. Liu, Y. Cheng, R. Li, and Z. Xu, “Investigation of the influence of the aberration induced by a plane interface on STED microscopy,” Opt. Express 17(3), 1714–1725 (2009). 18. X. Hao, C. Kuang, T. Wang, and X. Liu, “Effects of polarization on the de-excitation dark focal spot in STED microscope,” J. Opt. 12(115707) (2010). 19. Y. Xue, C. Kuang, S. Li, Z. Gu, and X. Liu, “Sharper fluorescent super-resolution spot generated by azimuthally polarized beam in STED microscopy,” Opt. Express 20(16), 17653–17666 (2012). 20. http://en.wikipedia.org/wiki/Circular polarization.


Introduction
Various types of far-field optical microscopes with superresolution have been proposed and realized for sub-diffraction-limit fluorescence imaging during the last decades, such as photoactive localization microscopy (PALM) [1], stochastic optical reconstruction microscopy (STORM) [2], structured illumination microscopy (SIM) [3], stimulated emission depletion (STED) microscopy [4], and others [5][6][7][8][9].In addition, mathematically calibrated microscopes such as fluorescence emission difference (FED) microscopy [10], scanning laser mode microscopy (SLAM) [11], and others [12] have attracted much interest recently for the ease with which they achieve superresolution.The essential key point of the latter methods is to subtract two fluorescence images captured by two excitation lasers with the same wavelength but separate solid and donut spots.However, the factors used in these microscopes to subtract the one from the other of two images are experiential.If the factor is too large, an exaggerated size reduction leading to even smaller than the real sample size may result.Meanwhile, some of the sample's important information may be missed due to over-processing.Therefore, a quantitative study of the parameters related to the subtraction factor must be conducted to obtain a definite criterion for determining appropriate subtraction factor.Previous reports [10][11][12] have mainly focused on the subtractions of the point spread function (PSF) without taking the sample structure into consideration.In fact, due to the removal of the negative values in the subtracted images, the data processing involved is nonlinear.As such, investigating only the PSF subtractions cannot sufficiently describe the real PSF of the subtraction microscope.In another case, due to the complexity of the convolution between the donut spot and different sample sizes, it is difficult to obtain an analytical description of the resolution in subtraction microscopes.
In this paper, a numerical investigation is carried out based on vector beam diffraction theory (VBDT) and two-dimensional convolution calculations.Using the calculated focal PSFs, the resolutions of the samples in the convoluted and subtracted images are analyzed and discussed along with the sample size, sample shape, numerical aperture (NA) of the microscope objective, and excitation laser field polarizations.Combined with the subtraction threshold line obtained from the simulation, the peak intensity after subtraction and resolution derivative versus the subtraction factor are found to be two useful parameters for evaluating the appropriate subtraction.

Principle and numerical models
Figure 1 shows the schematic principle of the fluorescence subtraction microscope.The sample is excited with the same laser wavelength but with two different focal spot shapes, which are usually a Gaussian distribution and a donut intensity distribution.The obtained fluorescence image by the latter is subtracted from that by the former as follows: I = FI G − γFI D , where FI G and FI D are the intensities of the two fluorescence images and γ is the subtraction factor.The produced image has a higher resolution compared with the conventional fluorescence image FI G .The most important point of this method is to obtain two suitable spot shapes and define a reasonable factor γ. Several combinations of spot shapes can be used for the subtraction microscope, such as the circular polarized Gaussian spot and vortex-phase-plate-generated donut spot [10,11], or the vector beam radially polarized solid spot and azimuthally polarized donut spot [12].The subtraction factor γ, which should be related to both the parameters of the two excitation lasers and the properties of the samples, plays a critical role in improving the resolution power.However, it can also introduce uncertainties.In the paper, an approximate numerical model is applied to investigate the image subtraction procedure and the factor γ. Because the sample is assumed to be very thin and excited without saturation, its structural shape along the axial direction can be disregarded.Thus, the intensity of the fluorescence is linearly proportional to the excitation laser intensity and to the density distribution of the fluorophores in the sample.The fluorescence excited by the two spots is assumed to have the same emission luminosity and collection efficiencies, and the detection PSFs and emission Stokes shifts of the fluorophores are not considered.Therefore, the fluorescence images can be described as the two-dimensional convolution of the excitation PSFs with the density distribution of the fluorophores in the sample.After the subtraction and the removal of negative values, the final intensity I(x, y) of the subtracted image can be expressed as follows: where I G and I D are the intensity-normalized PSFs of the excitation lasers, Ob is the shape (fluorophore density distribution) of the sample, γ is the subtraction factor, and * is the convo-lution.The PSFs are derived from the focal spots of the objective lens calculated by the VBDT [13][14][15][16][17][18][19][20] and expressed as follows: A L is the vector weight matrix of the lens, calculated as follows: The [P x , P y , P z ] is the polarization unit vector matrix, and has the values shown in Table 2.

Gauss Circ Rect
Ob(x, y) C 0 e The fluorophore density of the sample Ob is assumed to have a Gaussian distribution, and the densities with semicircular and rectangular types are also calculated for comparison.width at half maximum (FWHM), radius, and half width of the fluorophore density distributions ("Gauss", "Circ", and "Rect"), respectively.
The subtraction image I of the sample can be obtained by numerically calculating the PSFs and determining the convolution.The obtained FWHM(r, γ) of the sample in the subtraction images can be compared with the preset sample width r to evaluate the subtraction factor.In addition, the density types of the fluorohpore distributions, the NA of the objective lens, and the polarization combinations of the excitation lasers can also be investigated.Figure 2 shows the convolution and subtraction results with different sample sizes.Due to the convolution, the widths of the fluorescent distributions are broadened and the dark center of the donut spot is flattened to be non-zero compared with the initial PSFs, as shown in the 2D convoluted images in Fig. 2(a) and their central cut lines in Figs.2(c)-2(d).The sample is illustrated in Fig. 2(b).This illustrates that the donut center is sensitive to the sample sizes, and that the subtraction factor cannot be infinitely large because the FI G can be totally erased with some specific subtraction factors.This conclusion conflicts with the situation where only the PSFs are subtracted [10][11][12], in which case the subtraction factor γ can be infinitely large.Subtracting the PSFs and setting the negative values to be zero introduces additional nonlinearity and derivation, making the situation inappropriate.When subtracting with factors γ=0, 0.3, and 1.0, the obtained FWHM resolutions are 0.44λ , 0.33λ and 0.2λ with r=0.2λ ; and 0.49λ , 0.35λ and 0.13λ with r=0.3λ , respectively, as shown in Figs.2(e) and 2(f).This not only simply shows that a large subtraction factor leads to a better resolution, but also indicates that the FWHM decreases more rapidly as the subtraction factor increases in case of a large-radius sample.This can also be seen when the evolution of the FWHM resolution is compared with that of the subtraction factors in Fig. 2(i).The resolution with a sample size of r=0.3λ has a sharper slope than those with smaller sample sizes.Figures.2(e)-2(i) also show that with a smaller γ the sample may be resolved properly or below-subtraction.When the γ value is larger, the subtracted FWHM resolution of a large-size sample has an even smaller value than that of a small-size sample, which represents obvious over-subtraction.

Sample shape
The simulations with fluorophore densities of semicircular and rectangular types are carried out in comparison with the Gaussian type.The convoluted images and their cut lines using three sample shapes are shown in Figs.3(a), 3(b), and 3(d), 3(e), 3(h), 3(i), and the sample sizes are r=0.15λand 0.25λ , respectively.The samples with semicircular and rectangular shapes are illustrated in Fig. 3(c).They indicate that the contribution of three shapes to the FWHM resolution is considered approximately to be same with each other if the sample widths are small.When a sample size of r=0.15λ and factors of 0.3 and 1.0 are applied, the FWHM resolutions after subtraction with three different shapes are all close to 0.32λ and 0.21λ , respectively, as shown in Figs.3(f) and 3(g).This synchronized evolution of the FWHM resolution versus subtraction factors for three different types can similarly be found in Fig. 3(l).When the sample width increases to r=0.25λ , the FWHM curves of the three sample shapes are divergent, as shown in Figs.3(j), 3(k) and 3(m).The one with "Gauss" shape shows a smoother slope and it with the "Rect" shape shows a sharp decrease in resolution with the increase of the subtraction factors.Figures.3(n) and (o) give the FWHM resolution evolutions with three sample shapes versus the sample sizes with γ=0 and 0.4.They indicate that the FWHM resolution change with the subtraction factor is more sensitive to the sample sizes than to the sample shapes.For the sample sizes of r < 0.2λ in particular, the influence of the sample shapes on the resolution and subtraction factor can be ignored.

Subtraction factor
The principal criterion for evaluating the appropriate subtraction factor is that the obtained FWHM resolution in the subtracted images can be matched with the sample's real size.Figures.4(a)-4(c) give the 2D FWHM resolutions versus different sample sizes and subtraction factors with three sample shapes.In terms of the Gaussian-type samples, the obtained FWHM should be equal to the preset sample FWHM according to the subtraction criterion, which is displayed as the dotted black line in Fig. 4(a).The subtraction factors on this line are critical values for obtaining the right subtracted resolutions.The separate left and right regions are the below-and over-subtraction regions, respectively.This line can be picked out and plotted in Fig. 4(d) with the solid blue curve.We call it subtraction threshold line or γ 0 -r 0 line.With a specific subtraction factor γ 0 , the sample structures with sizes below the corresponding value r 0 can be imaged properly, and the sizes of r > r 0 are over-subtracted.The dotted blue threshold line in Fig. 4(d) is estimated simply by enhancing the FWHM with a factor of 1.4 for the confocal microscope case.However, using the unknown sample size to define the subtraction factor suitable for resolving the sample is not rational.Therefore, additional parameters combined  with the threshold line are needed to evaluate the proper subtraction.On the one hand, when Fig. 4(a) is coordinate transformed to Fig. 4(e), which is similar to Fig. 2(i), the FWHM slopes for the different sample sizes are clearly different and could be useful for identifying the real sample sizes.The absolute values of the slopes, defined as ΔFW HM/Δγ, are plotted with the threshold line (dotted black curve) in Fig. 4(f) with sample sizes of r=0.1λ , 0.2λ , 0.3λ , 0.4λ and 0.5λ .Although the threshold line is not systematic in the ΔFW HM/Δγ figure, relative sample sizes can still be compared among The larger the absolute value of ΔFW HM/Δγ, the larger the sample size is.On the other hand, the peak intensities of the sample in subtracted images decrease as the subtraction factor increases. Figure 4(g) shows the evolution map of the peak intensities, which are normalized to the initial peak intensities.The shape of the sample is also assumed to be Gaussian.The γ 0 -r 0 threshold line is plotted in Fig. 4(g).It is interesting to see that the γ 0 -r 0 line can horizontally and almost smoothly separate the below-and oversubtraction regions.The critical value associated with the peak intensity status for appropriate subtraction is around 0.6 even for different object sizes.This means that if the fluorescence peak intensity of the sample after subtraction drops below 60% of its initial value, the possibil-ity that it will be over-subtracted increases.It can be simply called I P factor.Although the best subtraction factor for the specific sample size varies, the decrease in peak intensity can be used to determine whether the factor is proper.The Fig. 4(h) gives an example to show that, when the γ=0.5 is used for r=0.3λ , the fluorescence peak intensity slightly drops below 0.6 and a preliminary over-subtraction appears.The corresponding fluorescence peak intensity and FWHM resolution after subtraction is 0.59 and 0.28λ , respectively.Similarly, the FWHM values after subtraction with semicircular and rectangular sample shapes should be approaches their preset FWHM sizes of √ 3r and 2r [illustrated in Fig. 3(c)].The threshold lines of them are shown in Figs.4(b) and 4(c) as the dotted black lines.Because these two types are rarely used in the imaging analysis, the results of them are not going to be shown in details.The resolution slopes of them are also identical for different sample sizes while their I P factors are found around 0.75 and 0.83, respectively.To inspect the effectiveness of the I P factor and the ΔFW HM/Δγ factor, further simulations are carried out with a laser wavelength of 488 nm, NA=1.4,and bead FWHMs of 200, 150, 100, and 50 nm with the Gaussian shape, respectively, shown in Fig. 5.After the 2D convolution with an intensity-normalized PSF, the FWHMs of them in the fluorescence image become 278, 244, 219, and 204 nm as shown in Fig. 5(a); the FI D and the sample is shown in Figs.5(b) and 5(c), respectively.With five different subtraction factors of γ=0.2, 0.4, 0.6, 0.8, and 1.0, the FWHM values of the four beads in the subtracted images decrease nonlinearly.Figure 5(d) displays the evolutions of their ΔFW HM/Δγ factors.That the sample with a larger absolute ΔFW HM/Δγ factor gets a larger size is in agreement with the prediction.Figure 5(e) shows the evolutions of their peak intensities after subtraction with the variation of the subtraction factor.The imaging with a sample width of 200 nm can be taken as an example.Its I P after subtraction drops below 0.6 when the subtraction factor goes further than 0.4.The FWHM resolution in the image reaches 163 nm as shown in Fig. 5(g) and this value decreases even further in Figs.5(h)-5(j) where they are all below the real width of this model sample.The resolutions of the other three samples are over-subtracted when their peak intensities drop below 60%.Therefore, the over-subtraction of the sample is predicted and demonstrated based on the peak intensity after subtraction.The over-subtracted areas are marked with dotted white circles in Figs.5(f)-5(j).As shown in Fig. 6, the simulations with an objective lens of NA=0.9 are also carried out using the same procedure.Because of the coarse resolution, the maps of the FWHM(r, γ) are upward shifted, and the threshold lines are also right shifted compared with Fig. 4(d).The value of the threshold line is also around 0.6.Therefore, the value of 0.6 could also be used as a criterion to inspect the over-subtraction.Figure 6(h) describes a critical subtraction with a Gaussian shape, sample radius of 0.3λ , and subtraction factor of 0.5.After subtraction, the peak intensity drops to 60% of the initial intensity, and the achieved resolution is around 0.3λ .

Polarization
The combination of radial and azimuthal polarizations is demonstrated to have fewer negative sidebands for subtraction microscopes [12].Simulations with this polarization combination and NA=1.4 are also realized here, as shown in Fig. 7.The FWHM resolution maps and ΔFW HM/Δγ factor are similar to those in Fig. 4 and are not discussed further.The I P factor is a little smaller than 0.6 in Fig. 7(h) but can still be considered flat.Figure 7(i) shows the right-subtraction status and implements the subtraction factor of 0.5 and sample FWHMs of 0.3λ .Although the polarization combination has changed, the criterion of the I P factor and ΔFW HM/Δγ still works.

Discussion
The simulations in this paper are carried out under several assumptions.First, the assumption of non-saturation allows the fluorescence obtained to be linearly proportional to the intensity of the excitation laser and the density of the fluorophores in the sample.In this case, the direct convolution of laser PSFs and the sample can be used to describe the conventional laser scanning fluorescence imaging.When a high excitation power is used, the fluorescence signals are intensity saturated, at which point the additional saturation factor should be included.Second, the detection PSF is assumed to be constant having the value of unity, which means that all of the fluorescence in one scanning step is collected.For confocal microscopy with different detection focal Airy disks and different fluorophores with different emission Stokes shifts, additional de-tection PSFs could be included to obtain more specific and accurate subtraction factors.Third, the samples are simulated with simple structures.For the biological samples with complex structures, the subtraction factors for the overlap regions still need to be addressed.In addition, dynamic subtraction factors could be a good measure for obtaining an ideal resolving power, but obtaining an accurate peak positioning requires the development of a complicated algorithm.However, the numerical simulation this paper can provide information of the FWHM resolutions as a function of the subtraction factor relevant to the sample size, shapes, NA, and polarizations in the subtraction microscope and to derive the appropriate subtraction threshold line.Meanwhile, this paper demonstrates that combined with the threshold line, the peak intensity after subtraction and the derivative of FWHM resolution are two useful parameters for evaluating the subtraction factor and estimating the real size of a sample.The samples with larger sizes are found to be sensitive to the subtraction factor and easily over-subtracted while the ones with smaller sizes can endure larger subtraction factor without over-processing.For the sample sizes of r < 0.2λ , the sample shapes do not play an important role in affecting the resolutions and determining the subtraction factors.When the fluorescence peak intensity in the subtracted image drops to 60% of the initial value, the corresponding subtraction factor is found to be properly chosen for the imaging even with different sample sizes.This criterion is also valid when a smaller NA is used, and can also roughly work when the radial and azimuthally polarized lights are used.The sample sizes can also be relatively compared according to the resolution derivative versus the subtraction factor.In summary, based on vector beam diffraction theory and two dimentional convolution, we have numerically invesigated the parameters that may affect the selection of subtraction factors in the subtraction microscopes and proposed that the fluorescence peak intensity after subtraction and resolution derivative can be possibly used as the key parameters for determing the appropriate subtraction factors and obtaining the sample sizes in the subtraction microscopes.

Table 1 .
AMP and A PHASE are the real amplitude and phase terms of the input beams, respectively.If the plane wave is used, then A AMP =1.Furthermore, A PHASE = e iϕ or A PHASE = 1 where a vortex phase plate is present or absent, respectively.Two groups of I G and I D are calculated according to the setting values listed in Table1.The second group is used only in Section 3.5.Setting values of I G and I D .

Table 3 .
Numerical models for the sample in Equation (1).